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Abstract 


We  have  assembled  a  data  set  of  nearly  1000  events  across  Eurasia  and  North  Africa  whose 
locations  and  origin  times  are  known  with  uncommon  accuracy.  These  “Ground  Truth”  (GT) 
events,  which  occur  in  23  event  clusters,  the  groomed  P  and  Pn  arrival  time  data  observed  at 
regional  distances  (A  <  20°),  and  about  1000  empirical  phase  path  anomalies  that  are 
constructed  from  each  of  the  clusters,  create  a  GT  data  base  to  be  used  to  assess  the  quality  of 
3-D  models  and  their  applicability  to  regional  locations.  We  assess  several  seismic  models  of  the 
crust  and  upper  mantle,  concentrating  on  the  1-D  model  AK135  and  the  3-D  vs  model  CUB2.0. 
Two  assessments  are  performed:  (1)  the  fit  of  model  predicted  travel  times  to  the  empirical  phase 
path  anomalies  and  (2)  the  ability  to  locate  the  GT  events.  The  3-D  vp  model  converted  from 
CUB2.0  vs  using  the  thermoelastic  properties  of  an  assumed  mantle  composition  fits  regional  P 
and  Pn  travel  times  remarkably  well,  both  in  geographical  pattern  and  absolute  level.  Rms-misfit 
is  about  1.0  sec  for  explosions,  and  1.25  sec  and  1.55  sec  for  earthquakes  with  location  accuracies 
of  about  5  km  and  10  km,  respectively.  We  conclude  on  the  basis  of  the  location  assessments  that 
the  intrinsic  regional  location  accuracy  is  ~  5  km  for  the  3-D  model  and.~  10km  for  the  1-D 
model.  As  the  number  of  stations  reduces,  both  1-D  and  3-D  model  locations  degrade  and 
variance  increases,  but  even  with  only  5  recording  stations  the  3-D  model  delivers  better  locations 
for  about  75%  of  events.  Our  findings  underscore  the  importance  of  GT  data  bases  in  assessing 
3-D  models,  and  show  that  global  3-D  vs  models  provide  a  sound  basis  on  which  to  build 
improvements  in  regional  location  capabilities. 
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Summary 


The  purpose  of  this  work  is  to  consider  the  feasibility  of  the  use  of  3-D  global  shear- wave 
models  of  the  crust  and  upper  mantle,  constructed  predominantly  with  surface  wave  dispersion 
information,  to  locate  small,  regional  seismic  events.  The  principal  interest  has  been  to  determine 
how  well  this  kind  of  model,  when  converted  to  compressional  velocity,  is  able  to  fit  regional  P 
and  Pn  travel  times  and  locate  small  events  using  only  regional  data  (distances  less  than  ~  20°). 
Global  shear  velocity  models  have  been  improving  rapidly  over  the  past  decade,  and  if  these 
models  are  able  to  predict  regional  P- wave  travel  times  well  they  would  provide  a  solid 
foundation  in  which  to  introduce  empirical  path  corrections. 

To  address  this  question,  we  developed  a  unique  validating  data  set  of  Ground  Truth  (GT) 
event  locations,  groomed  arrival  time  data,  and  empirical  path  anomalies  as  well  as  a  new  P- wave 
model  of  the  crust  and  uppermost  mantle  beneath  Eurasia,  the  Middle  East,  and  North  Africa. 
Validation  exercises  consisted  of  tests  of  the  ability  of  the  model  to  fit  the  regional  empirical  path 
anomalies  and  to  locate  the  GT  events  using  only  regional  data  from  the  groomed  arrival  time 
data  set.  The  principal  results  of  these  tests  are  described  in  section  3  of  this  document. 

Most  of  the  research  supported  by  this  contract  has  been  published  or  is  in-press,  much  in  the 
refereed  literature.  Tables  1  and  2  list  the  principal  publications  that  describe  this  research. 

The  work  supported  by  this  contract  divided  into  four  key  areas: 

(a)  model  construction, 

(b)  GT  data  base  development, 

(c)  ray  tracer  validation,  and 

(d)  model  validation. 

Tables  1  and  2  list  the  15  principal  publications  supported  directly  by  this  contract. 

We  constructed  two  models  of  the  crust  and  upper  mantle  beneath  Eurasia,  the  Middle  East,  and 
North  Africa  that  were  the  basis  for  our  tests.  Both  rest  on  new  group  velocity  measurements 
added  during  the  course  of  the  contract  as  well  as  phase  velocity  measurements  contributed  by 
Harvard  University  and  the  University  of  Utrecht.  The  models  differ  in  three  principal  respects: 
the  starting  crustal  model,  the  mode  of  surface  wave  tomography,  and  the  vs-to-vp  conversion 
procedure.  The  first  model,  called  CUB1.0,  is  built  on  CRUST5.1  (Mooney  et  al.,  1998),  uses 
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Table  1.  Publications  in  refereed  journals. 


Area 

Reference 

Model  Development 

tomography 

Barmin  et  al.  (2001) 

tomography 

Ritzwoller  et  al.  (2002b) 

inversion 

Shapiro  and  Ritzwoller  (2002a) 

inversion 

Villasenor  et  al.  (2001) 

data  assimilation 

Shapiro  and  Ritzwoller  (2002b) 

data  assimilation 

Engdahl  and  Ritzwoller  (2001) 

data  assimilation 

Ritzwoller  et  al.  (2002a) 

Ray  Tracer  Validation 

Villasenor  et  al  (2002) 

Model  Validation 

Levshin  and  Ritzwoller  (2002) 

Ritzwoller  et  al.  (2002c) 

Table  2.  Publications  in  unrefereed  journals  and  proceedings  volumes. 


Area  Reference 


Model  Development 

inversion 

Ritzwoller  et  al.  (1999) 

data  assimilation 

Ritzwoller  et  al  (2000) 

error  analysis 

Shapiro  et  al.  (2001) 

Ray  Tracer  Validation 

Villasenor  et  al  (2000) 

GT  Data  Base  Development 

Engdahl  et  al.  (2002) 

“Gaussian  tomography”  which  is  theoretically  based  on  geometrical  ray  theory  with  ad-hoc 
smoothing  constraints,  and  uses  an  empirical  logarithmic  scaling  relation  to  convert  from  vs  to  vp 
in  the  upper  mantle.  The  second  model,  called  CUB2.0,  is  based  on  CRUST2.0  (Laske,  personal 
communication),  “diffraction  tomography”  which  models  the  spatial  extent  of  surface  wave 
Fresnel  zones,  and  a  theoretical  conversion  between  vs  and  vp  based  on  recent  results  from 
mineral  physics.  We  show  below  that  CUB2.0  performs  better  in  the  validation  tests  than 
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CUB1.0.  We  also  investigated  the  addition  of  other  types  of  data,  including  Pn  and  Sn 
tomography,  heat  flow  data,  and  receiver  functions,  to  improve  the  model,  but  did  not  assimilate 
this  information  into  the  models  we  validated. 

We  tested  both  finite  difference  and  shooting  methods  to  trace  rays  and  compute  travel  times 
through  our  models.  The  codes  we  developed  agreed  quite  well  and  we,  ultimately,  chose  to  base 
our  validation  tests  on  the  shooting  method  because  it  is  much  faster  and  more  flexible  in  defining 
a  set  of  different  phases. 

The  assessment  of  a  seismic  velocity  model  requires  a  validating  data  set  that  is  independent  of 
the  data  used  to  derive  the  model.  For  this  study,  we  assembled  a  data  set  of  seismic  events  with 
magnitudes  above  the  IDC  event  definition  threshold  (currently  ~  3.75)  whose  locations  and 
origin  times  are  known  with  greater  than  usual  accuracy,  either  because  they  have  been  located 
with  a  local  seismic  network  or  because  they  were  anthropogenic  explosions.  These  so-called 
“Ground  Truth”  (GT)  events  are  used  to  obtain  reliable  empirical  estimates  of  source-station 
path  anomalies  (relative  to  a  1-D  reference  model)  that  can  be  compared  with  predictions  from 
seismic  models.  The  known  locations  of  the  GT  events  can  also  be  compared  to  locations 
obtained  from  seismic  models. 

The  GT  data  base  that  we  developed  consists  of  three  data  sets  that  are  useful  for  assessing  3-D 
seismic  models.  First,  there  is  the  set  of  989  GT  locations  that  are  available  to  test  location 
capabilities.  In  fact,  only  753  of  these  events  are  identified  as  GT5  or  better,  and  it  is  these  events 
that  proved  most  useful  in  this  study.  Second,  there  is  the  groomed  travel  time  data,  particularly 
regional  P  and  Pn,  from  the  ISC  and  NEIS  bulletins  as  well  as  some  data  from  regional  networks, 
that  formed  the  basis  for  the  HDC  analyses.  These  data  are  used  to  relocate  the  GT  events  in 
section  3  below.  Finally,  there  are  the  empirical  phase  path  anomalies  that  in  section  3  are 
compared  with  the  travel  times  predicted  from  several  models  (both  1-D  and  3-D). 

Studies  aimed  at  evaluating  the  location  capabilities  of  certain  models  or  location  techniques 
commonly  concentrate  exclusively  on  the  accuracy  of  the  locations  relative  to  some  benchmark. 
Although  this  is  also  an  important  part  of  the  present  study,  we  believe  that  the  ability  to  model 
regional  travel  time  data  is  a  better  measure  of  the  quality  of  3-D  models  and  is,  in  fact,  a  more 
robust  predictor  of  regional  location  capabilities.  This  is  because  evaluations  based  on  location 
alone  are  complicated  by  a  variety  of  factors  that  inhibit  clear  assessment;  e.g.,  variations  in 
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network  geometry  from  one  region  to  another,  different  mixes  of  regional  phases,  differential 
quality  of  reported  travel  times,  and  so  forth.  The  carefully  constructed  empirical  phase  path 
anomaly  data  set  gives  us  confidence  that  assessments  of  the  fit  to  regional  travel  time  data 
provide  meaningful,  relatively  unambiguous  information  about  the  quality  of  the  models. 

We  assessed  several  seismic  models  of  the  crust  and  upper  mantle,  concentrating  on  the  1-D 
model  AK135  and  the  3-D  vs  model  CUB2.0.  Two  assessments  are  performed: 

(1)  the  fit  of  model  predicted  travel  times  to  the  empirical  phase  path  anomalies  and 

(2)  the  ability  to  locate  the  GT  events. 

The  3-D  vp  model  converted  from  CUB2.0  vs  using  the  thermoelastic  properties  of  an  assumed 
mantle  composition  fits  regional  P  and  Pn  travel  times  remarkably  well,  both  in  geographical 
pattern  and  absolute  level.  Rms-misfit  is  about  1.0  sec  for  explosions,  and  1.25  sec  and  1.55  sec 
for  earthquakes  whose  locations  are  known  to  about  5  km  and  10  km,  respectively.  We  conclude 
on  the  basis  of  the  location  assessments  that  the  intrinsic  regional  location  accuracy  is  ~  5  km  for 
the  3-D  model  and  ~  10km  for  the  1-D  model.  As  the  number  of  stations  reduces,  both  1-D  and 
3-D  model  locations  degrade  and  variance  increases,  but  even  with  only  5  recording  stations  the 
3-D  model  deliver  better  locations  for  about  75%  of  events. 

Our  findings  underscore  the  importance  of  GT  data  bases  in  assessing  3-D  models,  and  show  that 
global  3-D  vs  models  provide  a  sound  basis  on  which  to  build  improvements  in  regional  location 
capabilities. 


1  Introduction 

The  seismic  location  inverse  problem  has  a  history  dating  back  almost  a  century,  at  least  to 
Geiger  (1912).  In  recent  years  there  has  been  a  resurgence  of  interest  in  improving  abilities  to 
estimate  reliable  seismic  locations  (earthquakes,  explosions)  (e.g.,  Thurber  &  Rabinowitz,  2000). 
There  are  several  reasons  for  this.  First,  there  is  the  desire  to  apply  advances  in  seismological 
models  and  practice  to  improve  the  locations  contained  in  global  catalogs  produced,  for  example, 
by  the  International  Seismological  Centre  (ISC),  the  U.S.  National  Earthquake  Information 
Service  (NEIS),  or  the  International  Data  Center  (IDC)  in  Vienna.  Second,  there  is  greater 
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emphasis  now  on  small  earthquake  location  that  has  been  generated  by  the  Comprehensive 
Nuclear  Test-Ban  Treaty  (CTBT),  which  specifies  a  location  uncertainty  of  1000  km2  or  less  for 
events  with  magnitudes  greater  than  or  equal  to  4.0  (e.g.,  Ringdal  &  Kennett,  2001).  Third,  there 
has  been  growing  recognition  of  the  importance  of  high  quality  locations  in  seismic  tomography 
and  fault  zone  characterization.  Efforts  to  improve  the  location  of  seismic  events  have  fallen  into 
two  broad  categories:  (1)  advances  in  travel  time  models  to  improve  absolute  locations  and  (2) 
the  development  of  multiple  event  location  methods  to  improve  relative  location  capabilities. 
These  categories  are  not  entirely  disjoint  because  if  a  cluster  of  relatively  located  events  contains 
one  or  more  reference  events  with  well-constrained  absolute  locations,  the  relatively  located 
events  can  be  tied  to  absolute  locations. 

Attempts  to  improve  body  wave  travel  time  models  divide  into  teleseismic  and  regional  studies, 
where  regional  distances  are  generally  defined  as  less  than  about  20°.  Recent  efforts  to  improve 
teleseismic  location  capabilities  have  concentrated  on  the  application  of  current-generation  1-D 
and  3-D  mantle  models.  For  example,  Engdahl  et  al  (1998)  showed  that  global  improvements  in 
teleseismic  locations  are  possible  by  applying  new  methods  using  the  1-D  model  AK135  (Kennet 
et  aZ.,  1995).  Several  researchers  have  investigated  how  the  use  of  3-D  tomographic  models  of  the 
whole  mantle  affects  hypocenter  determinations  based  on  teleseismic  travel  time  data  (e.g.,  Smith 
&  Ekstrom  (1996),  Antolik  et  al  (2001),  Chen  and  Willemann  (2001)).  These  studies  have 
established  that  teleseismic  mislocations  are  reduced  relative  to  1-D  model  locations  (e.g., 
Jeffreys-Bullen  travel-time  tables  (Jeffreys  &  Bullen,  1940);  PREM  (Dziewonski  &  Anderson, 
1981),  IASP91  (Kennett  &  Engdahl,  1991)  when  recent  generation  3-D  mantle  models  based  on 
either  global  (e.g.,  Su  et  a/.,  1994;  Su  &  Dziewonski,  1997)  or  local  basis  functions  (e.g.,  van  der 
Hilst  et  al. ,  1997;  Bijwaard  et  al ,  1998;  Boschi  et  al.  1999)  are  applied.  Efforts  to  predict  travel 
times  with  empirical  teleseismic  heterogeneity  corrections  have  been  more  rare  (e.g.,  Piromallo  & 
Morelli,  2001).  At  regional  distances,  travel  time  predictions  depend  on  the  ability  to  model  the 
crust  and  upper  mantle.  Regional  travel  time  predictions  are  less  accurate  than  teleseismic 
predictions  because  heterogeneity  is  much  stronger  and  because  data  that  constrain  the  crust  and 
upper  mantle  are  more  limited  in  distribution  and  somewhat  more  difficult  to  interpret 
unambiguously.  In  addition,  regional  triplication  crossover  distances  are  highly  variable  which 
makes  it  difficult  to  identify  phases  correctly.  This  is  true  both  for  Pn  triplications,  as  well  as  for 
distinguishing  Pn  from  P  arrivals.  For  these  reasons,  model-based  approaches  to  predicting 
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regional  travel  times  (e.g.,  Firbas,  2000)  have  played  a  subordinate  role  to  the  development  of 
purely  empirical  corrections  (e.g.,  Myers  &  Schultz,  2000)  that  are  combined  with  spatial 
interpolation  schemes  (e.g.,  Schultz  et  al. ,  1998). 

The  result,  however,  is  that  the  ability  of  3-D  models  to  predict  regional  travel  times  remains 
largely  uninvestigated  and,  consequently,  the  applicability  of  3-D  models  to  regional  location 
exercises  is  poorly  understood.  The  purpose  of  this  study  is  to  fill  in  this  gap  in  knowledge.  We 
are  able  to  take  this  step  for  two  principal  reasons.  First,  we  have  developed  a  global  3-D  model 
of  the  crust  and  uppermost  mantle  that  is  designed  to  be  used  to  predict  regional  travel  times. 

We  are  unaware  of  any  other  3-D  global  model  that  has  been  designed  for  this  purpose.  To 
predict  reliable  regional  travel  times  from  a  global  model  requires  careful  attention  to  the  crustal 
part  of  the  model  in  order  to  reduce  the  crust-mantle  trade-off  that  may  bias  the  mantle  part  of 
the  model.  In  addition,  regional  travel  times  are  controlled  largely  by  the  turning  point  of  the  ray. 
The  vertical  gradient  of  the  model,  therefore,  must  be  well  determined.  For  regional  travel  times, 
the  vertical  resolution  is  probably  more  important  than  the  lateral  resolution  of  the  model. 
Second,  we  have  also  developed  a  unique  data  base  of  nearly  1000  “Ground  Truth”  (GT) 
locations  in  23  event  clusters  across  much  of  Eurasia.  These  locations,  together  with  the  groomed 
arrival  time  data  set  and  the  empirical  phase  path  anomalies  that  emerge  for  each  cluster,  provide 
invaluable  information  needed  to  determine  the  capability  of  any  model  or  set  of  empirical  travel 
time  correction  surfaces  to  predict  regional  travel  times  and  locate  seismic  events  using  regional 
travel  time  data  alone. 

Improving  the  ability  to  locate  regional  events  is  important  for  two  principal  reasons.  First,  we 
believe  that  major  improvements  in  global  seismicity  catalogs  will  derive  principally  through 
better  use  of  regional  data,  at  least  for  continental  crustal  events.  As  seismic  networks  proliferate, 
the  quantity  of  regional  data  increases,  but  regional  travel  times  are  severely  down-weighted  in 
most  global  catalogs  due  to  the  inability  of  1-D  models  to  fit  these  phases  well  and  difficulties  to 
identify  phases  reliably.  The  need  for  consistency  in  global  catalogs  has  inhibited  the  adoption  of 
modern  models  largely  because  3-D  model  locations  using  teleseismic  data  provide  only  a 
moderate  improvement  over  locations  based  on  1-D  models.  Thus,  the  ISC  and  the  NEIS  continue 
to  use  the  Jeffreys-Bullen  travel-time  tables  and  the  IDC  uses  the  1-D  model  IASP91  (Kennett  h 
Engdahl,  1991).  Second,  the  concentration  on  monitoring  efforts  stimulated  by  the  CTBT  focuses 
attention  on  small  events  that  may  be  clandestine  nuclear  explosions.  Small  events  are  recorded 
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primarily  at  regional  distances,  perhaps  with  a  few  reported  teleseismic  phases.  Moreover,  the 
on-site  inspection  provisions  of  the  treaty  require  the  locations  of  these  small  events  to  be  of  high 
quality;  semi-major  axes  of  estimated  error  ellipses  are  specified  to  be  less  than  18  km.  The  focus 
of  research  efforts  motivated  by  the  technical  requirements  of  the  CTBT,  therefore,  is  on  regional 
locations  with  a  sparse  network  of  recording  stations. 

In  section  2,  we  discuss  the  procedure  used  to  construct  the  data  base  of  GT  locations,  the 
groomed  arrival  time  data  set,  and  empirical  phase  path  anomalies  for  each  event  cluster.  We  also 
describe  there  the  grid-search  location  procedure  to  locate  the  GT  events  using  regional  data 
alone.  The  3-D  model  construction  and  the  ray-tracing  procedure  are  discussed  in  section  3.2. 

The  first  test  of  the  3-D  model  is  its  ability  to  predict  regional  travel  times.  We  estimate  in 
section  3.3  the  fit  to  the  empirical  phase  path  anomalies.  The  second  test  of  the  3-D  model  is  the 
ability  to  locate  the  GT  events  using  regional  data  alone.  We  present  in  section  3.4  results  of  the 
grid-search  location  for  21  event  clusters.  We  restrict  the  study  to  the  use  of  P  and  Pn  travel 
times  for  stations  within  20°  of  the  source  locations.  While  optimal  regional  locations  would 
include  crustal  phases,  later  P  phases,  and  S  phases,  the  consideration  of  these  phases  is  not 
necessary  to  test  the  performance  of  the  3-D  model  and,  in  fact,  would  complicate  the 
interpretation  of  the  location  exercise. 

2  Methods  and  Procedures 

2.1  Motivation 

The  assessment  of  a  seismic  velocity  model  requires  a  validating  data  set  that  is  independent  of 
the  data  used  to  derive  the  model.  For  this  study,  we  have  assembled  a  data  set  of  seismic  events 
with  magnitudes  above  the  IDC  event  definition  threshold  (currently  m &  ~  3.75)  whose  locations 
and  origin  times  are  known  with  greater  than  usual  accuracy,  either  because  they  have  been 
located  with  a  local  seismic  network  or  because  they  were  anthropogenic  explosions.  These 
so-called  “Ground  Truth”  (GT)  events  are  used  to  obtain  reliable  empirical  estimates  of 
source-station  path  anomalies  (relative  to  a  1-D  reference  model)  that  can  be  compared  with 
predictions  from  seismic  models.  The  known  locations  of  the  GT  events  can  also  be  compared  to 
locations  obtained  from  seismic  models. 
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Figure  1:  Locations  of  the  explosion  and  earthquake  clusters.  Numbers  refer  to  clusters  identified 
in  Table  3.  Lines  A  -  A1  and  B  —  B*  refer  to  the  vertical  model  profiles  shown  in  Figure  9. 

In  this  section,  we  discuss  the  development  of  the  validating  data  set  by  multiple  event  relocation 
of  clustered  earthquake  and  explosion  sequences  in  Eurasia  and  northern  Africa.  Primarily,  we 
use  phase  arrival  time  data  reported  to  the  ISC  and  to  the  NEIS,  but  for  some  clusters  we  have 
obtained  additional  data  (usually  at  regional  distances)  from  other  sources.  The  arrival  time  data 
are  first  “groomed”  on  an  event  by  event  basis,  using  a  procedure  described  by  Engdahl  et  al 
(1998)  in  the  development  of  the  EHB  catalog.  A  multiple  event  location  method  is  then  used  to 
refine  the  locations  and  identify  out-liers  in  the  arrival  time  data.  The  resulting  catalog  of  GT 
events,  the  groomed  travel  time  data  base,  and  the  empirical  phase  path  anomalies  create  a  data 
base  that  can  be  used  in  experiments  to  assess  3-D  models  and  to  address  systematic  errors  in 
event  location. 

The  Ground  Truth  events  occur  in  the  clusters  shown  in  Figure  1  and  are  listed  in  Table  3. 
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Table  3.  Cluster  parameters  for  three  categories  of  events:  (1)  explosions,  (2)  GT5 

earthquakes,  (3)  GT10  earthquakes. 


Cluster  No. 

Lata 

LonQ 

Depth0 

A  t\ 

nevtc 

nref* 

GTze 

Name 

1 

47.875 

48.139 

0.1 

0.53 

7 

7 

GT1 

Azgir 

2 

49.954 

78.871 

0.0 

0.81 

100 

100 

GT1-2 

Balapan 

3 

49.784 

78.072 

0.0 

0.66 

151 

146 

GT1-2 

Degelen 

4 

41.580 

88.605 

0.0 

0.43 

20 

13 

GT1-2 

Lop  Nor 

5 

73.350 

54.820 

0.0 

-1.65 

29 

28 

GT1-2 

Novaya 

6 

24.049 

5.040 

0.0 

-0.08 

5 

5 

GT0 

Sahara 

7 

36.938 

35.825 

33.6 

-2.89 

24 

4 

GTS 

Adana 

8 

23.486 

70.265 

19.6 

0.36 

107 

6 

GT5 

Bhuj 

9 

30.590 

79.124 

12.1 

-0.18 

86 

8 

GT5 

Chamoli 

10 

40.798 

31.219 

9.1 

-0.64 

41 

2 

GT5 

Duzce 

11 

39.588 

39.805 

6.4 

-1.73 

9 

3 

GT5 

Erzin 

12 

38.822 

70.560 

7.5 

-1.55 

28 

4 

GTS 

Garm 

13 

35.235 

-3.930 

6.2 

-0.91 

38 

3 

GT5 

Hoceima 

14 

40.740 

30.223 

10.4 

-0.51 

34 

5 

GT5 

Izmit 

15 

17.236 

73.748 

8.2 

-0.29 

31 

10 

GT5 

Koyna 

16 

42.475 

43.737 

7.0 

-1.28 

35 

5 

GTS 

Racha 

17 

57.092 

122.276 

31.5 

-1.53 

8 

3 

GTS 

Siberia 

18 

13.513 

51.069 

10.0 

-0.45 

55 

5 

GT10 

Aden 

19 

28.754 

34.608 

14.0 

-1.94 

36 

1 

GT10 

Aqaba 

20 

27.482 

33.864 

12.8 

-1.49 

27 

1 

GT10 

Gubal 

21 

39.660 

76.988 

19.6 

-2.75 

72 

1 

GT10 

Jiashi 

22 

40.957 

44.215 

5.6 

-1.60 

11 

2 

GT10 

Spitak 

23 

33.479 

57.111 

14.6 

-1.95 

35 

2 

GT10 

Tab  as 

Explosions 

312 

299 

GT5  Earthquakes 

441 

53 

GT10  Earthquakes 

236 

12 

Total 

989 

314 

a  Hypocentroid  location,  depth  in  km. 

b  Cluster  time  baseline  shift  applied  to  teleseismic  and  regional  data,  in  sec. 
c  Number  of  GT  events  in  the  cluster. 
d  Number  of  reference  events  in  the  cluster. 
e  95%  confidence  in  location  accuracy  better  than  x  km. 
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Figure  2:  Schematic  of  the  HDC  analysis  used  to  generate  the  GT  data  base  in  this  study,  (a) 
Relative  locations  of  cluster  event  epicenters  (“cluster  vectors”)  are  shown  with  90%  confidence 
ellipses,  (b)  The  hypocentroid  (solid  star)  is  located  using  all  cluster  event  arrival  time  data 
(appropriately  shifted).  The  box  corresponds  to  (a),  (c)  The  relative  locations  from  (a)  are  added 
to  the  absolute  location  of  the  hypocentroid  in  (b).  (d)  The  location  bias  is  removed  by  shifting 
the  cluster  to  match  the  reference  event  locations  .  Reference  locations  for  event  1  and  3  are  shown 
by  open  stars. 
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travel  time  model  to  locate  the  hypocentroid,  just  as  single  event  locations  are  biased.  The 
direction  and  degree  of  bias  depend  on  the  nature  of  departures  of  the  1-D  model  from  the  true 
velocity  structure,  combined  with  the  particular  distribution  of  observing  stations.  To  minimize 
this  bias,  we  shift  the  hypocentroid  (i.e.,  the  entire  cluster)  in  space  and  time,  as  illustrated  in 
Figure  2d,  in  order  to  provide  an  optimal  match  between  the  reference  event  locations  and  the 
corresponding  cluster  events.  This  brings  all  events  in  the  cluster  into  close  alignment  with 
ground  truth.  We  refer  to  the  resulting  set  of  989  locations  listed  in  Table  3  as  the  GT  locations. 
The  accuracy  of  the  GT  locations  and  the  shift  in  the  origin  time  (“cluster  time  baseline  shift”) 
are  discussed  further  below. 

For  several  reasons,  Pn  and  P  arrival  time  data  were  the  primary  phases  used  in  the  HDC 
analyses.  First,  they  are  ordinarily  read  as  first  arriving  phases  and  are  therefore  more  accurately 
picked.  Second,  at  regional  distances  secondary  phases  usually  arrive  in  the  Pn/ P  coda,  making 
them  not  only  difficult  to  pick  but  also  making  it  difficult  to  identify  the  phase  type  correctly. 
These  difficulties  arise  mainly  because,  for  most  of  the  clusters  studied,  the  propagation  paths  at 
regional  distances  for  these  phases  pass  through  regions  of  strong  lateral  heterogeneity,  leading  to 
complex  waveforms.  However,  there  are  exceptions,  such  as  propagation  across  the  Indian  Craton 
where  Sn  is  observed  as  a  coherent  phase  that  can  be  read  accurately;  Sn  was  used  in  the  analyses 
of  the  Bhuj  and  Chamoli  clusters. 

The  reference  event  source  parameters,  which  form  the  heart  of  the  GT  data  base,  were  compiled 
from  many  sources,  but  were  based  exclusively  on  work  done  by  the  original  investigators  in  each 
source  region.  We  attempt  to  validate  the  claimed  accuracy  and  assign  a  GT  level  of  accuracy  by 
analyzing  the  local  network  data  and  through  tests  of  consistency  with  the  HDC  results.  The 
epicenter  accuracy  of  candidate  reference  events  is  first  assessed  on  the  basis  of  station  geometry, 
described  by  fixed  confidence  bounds  using  a  Ground  Truth  GTXC%  classification,  where  the  “X” 
suffix  designates  location  accuracy  in  km  and  C%  is  the  percent  confidence  (Bondar  et  al. ,  2002). 
For  example,  events  that  are  accurate  to  within  5  km  at  a  95%  confidence  level  are  designated 
GT595%.  Bondar  et  al.  (2002)  argue  that  crustal  events  are  located  with  5  km  accuracy  or  better 
at  the  95%  confidence  level  if  they  are  located  with  (1)  at  least  10  stations  within  250  km,  (2)  an 
azimuthal  gap  less  than  110°,  (3)  a  secondary  azimuthal  gap  less  than  160°,  and  (4)  at  least  one 
station  within  30  km  of  the  epicenter.  The  final  constraint  provides  some  confidence  in  depth. 
Secondary  azimuthal  gap  is  defined  as  the  largest  azimuthal  gap  filled  by  a  single  station.  In 
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addition,  it  is  desirable  for  local  network  solutions  to  be  based  on  both  P—  and  S— wave  arrival 
times  and  a  local  velocity  model.  For  most  clusters  presented  in  this  study,  original 
documentation  that  included  station  locations,  velocity  model  used,  arrival  times,  and  details  of 
the  solutions  was  compiled  and  carefully  reviewed  according  to  the  above  criteria. 

We  use  “GT-level”  classifications  in  two  senses.  First,  is  the  sense  of  the  accuracy  level  of 
individual  reference  events,  as  discussed  above.  Second,  we  use  the  same  concept  to  describe  the 
level  of  accuracy  of  entire  clusters  after  shifting  to  best  match  the  reference  locations.  The 
GT-level  assigned  to  a  cluster  is  more  subjective,  related  to  the  GT-level  of  individual  reference 
events,  the  number  of  reference  events,  the  degree  of  consistency  between  the  relative  locations  as 
expressed  by  the  HDC  analysis  and  reference  data,  and  other  factors.  For  example,  the  Aqaba 
cluster  is  calibrated  by  a  single  reference  event,  which  qualifies  as  GT5,  but  we  assign  the  cluster 
a  GT10  level  of  accuracy  because  there  is  only  a  single  reference  event  and  because  the  event 
occurred  in  a  region  of  very  strong  lateral  heterogeneity,  the  transition  from  continental  to 
oceanic  lithosphere. 

During  the  course  of  this  study  we  discovered  that  many  proposed  reference  events  are  in  fact 
seriously  mislocated.  The  consistency  between  the  relative  locations  as  determined  by  HDC 
analysis  of  global  arrival  time  data  and  the  relative  locations  derived  from  the  reference  event 
data  is  one  of  the  criteria  we  use  to  validate  candidate  reference  events.  Shifts  in  epicenter  and 
origin  time  to  match  the  reference  event  locations  are  typically  in  the  range  of  5-15  km  and  ±2 
seconds,  respectively. 

We  shift  the  hypocentroid  in  space  and  time  to  achieve  the  best  fit  for  absolute  locations  and 
origin  times  between  the  HDC-derived  locations  and  the  reference  event  locations,  as  illustrated  in 
Figure  2d.  The  “cluster  time  baseline  shift7’  is  the  average  difference  between  the  reference  event 
origin  times  and  the  origin  times  obtained  by  HDC  analysis  for  those  same  events  (Table  3).  This 
correction  is  added  to  the  HDC-derived  origin  times  to  produce  a  best  estimate  of  absolute  origin 
times.  Positive  corrections  (i.e.,  shorter  travel  times)  indicate  faster  velocities  relative  to  AK135. 
The  cluster  time  baseline  shift  arises  in  part  from  differences  between  the  1-D  model  travel  time 
used  in  the  HDC  analysis  and  the  true  travel  times,  but  may  also  reflect  bias  in  the  reference 
event  origin  times.  Myers  et  al  (2001)  have  shown  that  origin  times  estimated  from  local  network 
data  can  be  biased  by  several  seconds  if  the  local  velocity  model  is  not  well-calibrated,  even 
though  the  hypocenters  may  be  estimated  with  high  accuracy.  The  cluster  time  baseline  shifts 
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assist  in  phase  identification  and  are  also  used  to  compare  results  between  clusters  and  are 
included  in  estimates  of  empirical  phase  path  anomalies  relative  to  AK135,  as  discussed  below. 

Our  HDC  analyses  are  less  afflicted  by  phase  misidentification  than  single  event  locations  for  two 
reasons.  First,  the  relative  locations  of  cluster  events  are  determined  by  differences  in  the  travel 
time  residuals,  not  their  absolute  value.  Second,  we  use  the  information  on  cluster  residuals  (after 
average  path  error  is  removed  and  reading  errors  have  been  estimated),  to  identify  out-liers  which 
may  be  caused  by  phase  identification  errors.  Phase  misidentification  remains  as  a  problem  for 
the  estimation  of  the  hypocentroid,  but  its  influence  is  reduced  by  the  much  larger  number  of 
data  used  (several  thousand),  compared  to  most  single  event  locations.  In  any  case,  the 
hypocentroid  is  ultimately  shifted  to  match  the  reference  event  locations. 

This  method  has  resulted  in  a  data  base  of  23  event  clusters,  including  6  explosion  clusters  with 
source  locations  generally  known  to  better  than  2  km,  and  17  earthquake  clusters,  11  of  whose 
locations  are  believed  to  be  accurate  to  within  5  km  and  the  remainder  to  within  10  km.  There 
are  989  GT  events  that  compose  these  23  clusters  of  which  753  are  GT5  or  better.  The  locations 
of  the  clusters  studied  are  plotted  in  Figure  1.  Relevant  parameters  for  all  clusters  are  listed  in 
Table  3.  The  location  and  depth  given  are  for  the  shifted  hypocentroid,  representing  the  best  fit 
to  reference  event  locations. 

2.3  The  Groomed  Arrival  Time  Data  Set 

A  significant  advantage  of  multiple  event  location,  relative  to  single  event  location,  is  the 
opportunity  to  perform  an  iterative  “grooming”  process  in  order  to  identify  and  remove  out-liers 
in  the  residuals  until  the  normalized  (by  current  assumed  reading  error),  de-meaned  residuals 
resemble  a  standard  normal  distribution.  This  results  in  a  data  set  of  groomed  arrival  times  for 
each  of  the  events  in  the  GT  data  base.  The  number  of  these  data  is  listed  in  Table  4  and  Figure 
3  shows  a  plot  of  the  Pn  and  P  residuals  for  all  groomed  arrival  times  included  in  this  data  base, 
adjusted  by  the  cluster  time  baseline  shift.  Not  unexpectedly,  estimates  of  the  residual  median 
and  spread  show  considerable  variation  over  the  regional  distance  range.  In  particular,  residuals 
for  Pn  arrivals  between  about  9  and  17  degrees,  corresponding  to  ray  paths  bottoming  in  the 
lithosphere  well  below  the  crust  and  perhaps  encountering  low  velocity  zones,  are  more  spread, 
with  a  maximum  spread  at  distances  of  12-13  degrees.  A  contributing  factor  to  this  spread  at 
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larger  regional  distances  is  the  difficulty  of  phase  identification;  e.g.,  distinguishing  between 
various  Pn  branches  or  between  Pn  and  P  in  regions  where  the  crossover  distance  between  the 
first-arriving  branches  of  these  phases  is  poorly  known. 


Table  4.  Number  of  travel  times  in  the  groomed  arrival  time  and  the  empirical  phase 
path  anomaly  data  sets. 


Groomed  Arrival  Times 

Empirical  Path  Anomalies 

Pn 

P 

Pn 

P 

Explosions 

1521 

477 

121 

35 

GT5 

8974 

1380 

558 

108 

GT10 

3106 

539 

157 

35 

total 

13601 

2396 

836 

178 

2.4  The  Empirical  Phase  Path  Anomalies 

The  basic  premise  of  all  multiple  event  location  methods  is  that  path  anomalies  from  each  station 
to  all  observed  events  in  a  given  cluster  are  correlated.  Thus,  multiple  event  location  analysis 
produces  robust  estimates  of  source-station  path  anomalies  that  are  far  more  difficult  to  extract 
from  single  event  location  catalogs.  The  set  of  groomed  residuals  for  regional  P  and  Pn,  relative 
to  shifted  hypocenters  derived  by  HDC  analysis  and  adjusted  for  cluster  time  baseline  shifts,  are 
used  to  calculate  source-station  phase  path  anomalies.  These  anomalies  are  estimated  relative  to 
the  1-D  reference  model  AK135.  Medians  and  spreads  for  all  groomed  residuals  for  each  cluster 
are  calculated  for  all  the  phases  of  interest  (Pg,  Pn,  P,  Sg,  Sn  and  S)  at  each  of  the  reporting 
stations  for  that  cluster,  from  single  event  location  catalogs.  The  set  of  groomed  residuals  for 
regional  P  and  Pn,  relative  to  shifted  hypocenters  derived  by  HDC  analysis  and  adjusted  for 
cluster  time  baseline  shifts,  are  used  to  calculate  source-station  phase  path  anomalies.  These 
anomalies  are  estimated  relative  to  the  1-D  reference  model  AK135.  Medians  and  spreads  for  all 
groomed  residuals  for  each  cluster  are  calculated  for  all  the  phases  of  interest  ( Pg ,  Pn,  P,  Sg.  Sn 
and  S)  at  each  of  the  reporting  stations  for  that  cluster.  The  resulting  source-station  “empirical 
phase  path  anomalies”  (the  median)  are  accepted  with  a  minimum  requirement  of  five 
observations  and  a  spread  of  less  than  1.40  sec  and  2.8  sec  for  P—  and  S— type  phases, 
respectively.  (“Spread”  is  a  robust  analog  to  the  standard  deviation.) 
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Figure  3:  Residual  density  plot  for  Pn  and  P  phases  in  the  groomed  arrival  time  data  set.  Residuals 
are  taken  with  respect  to  AK135  travel  times.  Station  cluster  residuals  corrected  for  the  cluster 
event  static  are  binned  at  0.5  degrees  in  distance  and  0.5  sec  in  time.  Bin  hit  counts  are  plotted 
according  to  the  scale.  Estimates  of  the  residual  median  and  spread  at  one-degree  intervals  are  also 
plotted  as  circles  and  vertical  lines,  respectively. 
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Figure  4:  Reduced  Pn  (white)  and  P  (black)  empirical  phase  path  anomalies  plotted  as  the  median 
and  spread  estimates  with  respect  to  distance.  Travel  times  are  corrected  for  cluster  time  baseline 
shifts,  (a)  Results  are  presented  as  residuals  with  respect  to  the  AK135  travel  times,  (b)  Results 
are  presented  near  the  Pn/P  cross-over  as  reduced  travel  times  (relative  to  8.0  km/s).  AK135 
predicted  Pn  and  P  travel  time  branches  are  also  shown. 

Ordinarily,  there  are  too  few  Pg  phases  in  the  data  base  and  the  5- type  phase  residuals  are  too 
noisy  to  obtain  useful  results.  At  regional  distances,  however,  there  are  836  Pn  and  178  P  path 
anomalies,  respectively.  Results  are  plotted  in  Figure  4.  These  empirical  estimates  range  from 
about  -7.5  sec  to  T5.0  sec,  with  the  largest  anomalies  occurring  at  distances  from  11  to  18  degrees 
(similar  to  Figure  3).  In  this  distance  range,  Pn  typically  has  a  triplication  and  P  has  a  back 
branch  (Figure  4b)  making  phase  identification  difficult. 

2.5  Location  Method 

We  developed  a  grid-search  location  method  to  be  used  in  the  second  assessment  of  the  3-D 
model.  In  addition  to  the  three  principal  unknowns,  origin  time  (io),  depth  (z),  and  epicentral 
location  (z,y),  we  consider  phase  identification  to  be  unknown.  We  will  fix  the  depth  because  it 
trades  off  with  origin  time  and  fixing  it  has  little  effect  on  the  epicentral  location.  To  simplify  the 
inversion,  we  re-identify  phases  only  at  the  grid-center  rather  than  at  each  grid  node  separately, 
and  use  only  Pn  and  P  data  observed  between  ~  3°  and  20°.  The  use  of  Sn  and  S  travel  times 
would  improve  the  locations,  particularly  for  locations  with  few  reporting  stations.  The  quality 
and  distribution  of  reported  regional  S- phases,  however,  are  much  more  variable  than  P-phases 
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and  their  use  would  make  the  results  more  difficult  to  interpret.  These  simplifications  result  in 
much  more  stable  location  results,  but  the  location  method  differs  substantially  from  procedures 
that  are  used  to  construct  global  catalogs.  Thus,  our  reports  of  location  capabilities  are  probably 
more  meaningful  in  a  relative  sense  (e.g.,  location  error  with  the  3-D  model  relative  to  the  1-D 
model  location)  than  an  absolute  sense.  The  3-D  model  that  is  tested  here  is  CUB2.CLTH  and  the 
1-D  model  is  AK135. 

The  location  procedure  progresses  in  five  steps.  (1)  Fix  event  depth.  The  depth  z  is  fixed  to 
the  value  in  the  GT  location  data  base.  (2)  Create  the  grid.  The  grid  is  chosen  with  the 
grid-center  at  the  EHB  location  for  earthquakes  and  the  PDE  location  for  explosions,  and  nodes 
are  located  every  1  km  to  create  a  50  km  x  50  km  grid.  The  PDE  location  is  chosen  for  the 
grid-center  for  explosions  because  the  EHB  location  for  explosions  is  commonly  of  GTO-2  quality, 
which  would  create  an  asymmetry  between  earthquakes  and  explosions.  The  EHB  and  PDE 
locations  are  probably  more  accurately  characterized  as  GT10  -  GT15.  (3)  Phase  identification 
and  origin  time.  With  a  hypothesized  epicenter  at  the  grid-center,  each  travel  time  is  identified 
as  either  P  or  Pn  and  the  origin  time  is  shifted  (relative  to  an  input  estimate)  to  give  the  smallest 
over-all  rms-residual.  (4)  Out-lier  rejection.  Residuals  larger  than  ±3  sec  using  the  3-D  model 
are  considered  to  be  out-liers  and  are  rejected.  The  3  sec  criterion  is  based  on  fits  to  the  empirical 
phase  path  anomalies  (e.g.,  Figure  17,  discussed  later  in  the  report). 

It  is  about  the  3a  level  for  explosions  and  2cr  for  earthquakes.  The  details  of  this  choice  have  little 
effect  on  the  overall  statistics  because  the  number  of  observations  that  are  rejected  is  small.  In 
principal,  this  procedure  could  be  performed  separately  at  each  grid  node.  This  would,  however, 
produce  different  data  sets  for  each  model,  which  would  complicate  interpretation  of  the  results. 
We  also  believe  that  the  3-D  model  accurately  identifies  erroneous  measurements.  (5)  Epicenter 
estimate.  Hypothetical  epicenters  are  moved  to  each  grid  node  and  the  estimated  event  location 
is  identified  with  the  grid  node  that  produces  the  minimum  misfit  to  the  observed  travel  times  in  ' 
the  groomed  arrival  time  data  set. 

An  example  of  the  location  grid  for  an  explosion  on  the  Lop  Nor  test  site  is  shown  in  Figure  5. 
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Figure  5:  Example  of  the  grid-search  location  algorithm.  Contours  or  rms-misfit  using  a  3-D  model 
for  a  grid  of  hypothesized  epicenters  for  an  explosion  at  the  Lop  Nor  test  site  (May  26,  1990;  54 
reporting  stations).  Units  are  seconds.  The  GT1-2  location  is  indicated  with  a  star,  the  best-fit 
locations  using  the  3-D  and  1-D  models  are  shown  with  a  triangle  and  a  circle,  respectively,  and 
the  PDE  location  is  shown  with  a  diamond. 
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3  Results  and  Discussion 


3.1  Evaluation  Metric 

The  method  discussed  in  section  2  has  produced  three  data  sets  that  are  useful  for  assessing  3-D 
seismic  models.  First,  there  is  the  set  of  989  GT  locations  that  are  available  to  test  location 
capabilities.  In  fact,  only  753  of  these  events  are  identified  as  GT5  or  better,  and  it  is  these  events 
that  will  be  most  useful  in  this  study.  Second,  there  is  the  groomed  travel  time  data,  particularly 
regional  P  and  P„,  from  the  ISC  and  NEIS  bulletins  as  well  as  some  data  from  regional  networks, 
that  formed  the  basis  for  the  HDC  analyses.  These  data  will  be  used  to  relocate  the  GT  events  in 
section  3.4  below.  Finally,  there  are  the  empirical  phase  path  anomalies  that  in  section  3.3  are 
compared  with  the  travel  times  predicted  from  several  models  (both  1-D  and  3-D). 

Studies  aimed  at  evaluating  the  location  capabilities  of  certain  models  or  location  techniques 
commonly  concentrate  exclusively  on  the  accuracy  of  the  locations  relative  to  some  benchmark. 
Although  this  is  also  an  important  part  of  the  present  study,  we  argue  that  the  ability  to  model 
regional  travel  time  data  is  a  better  measure  of  the  quality  of  3-D  models  and  is,  in  fact,  a  more 
robust  predictor  of  regional  location  capabilities.  This  is  because  evaluations  based  on  location 
alone  are  complicated  by  a  variety  of  factors  that  inhibit  clear  assessment;  e.g.,  variations  in 
network  geometry  from  one  region  to  another,  different  mixes  of  regional  phases,  differential 
quality  of  reported  travel  times,  and  so  forth.  The  carefully  constructed  empirical  phase  path 
anomaly  data  set  gives  us  confidence  that  assessments  of  the  fit  to  regional  travel  time  data 
provide  meaningful,  relatively  unambiguous  information  about  the  quality  of  the  models.  The 
groomed  travel  time  data  set  alone  is  probably  too  noisy  to  be  used  for  this  purpose. 

3.2  3D  Model  and  Travel  Time  Computation 
3.2.1  Model  Construction 

The  3-D  model  is  based  on  broad-band  surface  wave  group  and  phase  speed  measurements.  The 
group  velocities  were  measured  using  the  method  described  by  Ritzwoller  h  Levshin  (1998), 
which  is  a  frequency  -  time  method  that  involves  analyst  interaction  to  choose  the  frequency  band 
of  measurement  and  to  guide  the  extraction  of  the  fundamental  mode  from  noise,  scattered  and 
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multipathed  signals,  overtones,  and  fundamental  modes  of  different  wave  types.  We  use  group 
speed  measurements  from  16  s  to  200  s  period  for  Rayleigh  waves  and  from  16  s  to  150  s  period 
for  Love  waves.  The  phase  speed  measurements  were  performed  at  Harvard  University  (Ekstrom 
et  al,  1997)  and  Utrecht  University  (Trampert  &  Woodhouse,  1995)  separately  and  we  merged 
these  data  sets.  The  phase  speed  measurements  extend  from  40  s  to  150  s  for  both  Rayleigh  and 
Love  waves.  We  use  measurements  only  from  earthquakes  shallower  than  50  km  to  reduce  the  size 
of  the  source  group  time  shifts,  which  we  do  not  attempt  to  correct  (Levshin  et  al.,  1999).  All 
measurements  are  subjected  to  the  quality  control  procedures  described  by  Ritzwoller  &  Levshin 
(1998),  but  the  number  of  group  speed  measurements  has  multiplied  several  times  since  that 
study. 

The  construction  of  the  group  and  phase  speed  maps  is  founded  on  the  tomographic  method  of 
Barmin  et  al.  (2001),  which  is  based  on  geometrical  ray-theory  with  intuitive  Gaussian  smoothing 
constraints  to  simulate  surface  wave  sensitivities.  We  refer  to  this  method  as  Gaussian 
tomography.  In  fact  we  apply  an  update  of  this  method  (Ritzwoller  et  al.,  2002),  referred  to  as 
diffraction  tomography,  which  uses  a  simplified  version  of  the  scattering  sensitivity  kernels  that 
emerge  from  the  Born  or  Rytov  approximations.  Diffraction  tomography  accounts  for 
path-length-dependent  sensitivity,  wave-front  healing  and  associated  diffraction  effects,  and 
provides  a  more  accurate  assessment  of  spatially  variable  resolution  than  traditional  tomographic 
methods.  The  resolution  procedure  is  described  by  Barmin  et  al.  (2001).  We  produce  a  resolution 
surface  at  every  nodal  point  on  the  globe,  fit  a  cone  to  the  surface  in  the  neighborhood  of  each 
node,  and  define  resolution  as  the  half-width  of  the  base  of  the  cone  (or  identically,  the  full-width 
at  half- max).  Surface  wave  resolution  estimates  averaged  over  the  region  of  study  are  presented  in 
Figure  6.  Surface  wave  Fresnel  zones  widen  appreciably  at  long  periods  and  narrow  only  near 
sources  and  receivers,  so  resolution  is  best  at  short  periods  and  in  areas  with  sources  or  receivers. 

The  dispersion  maps  are  produced  on  a  2°  x  2°  grid  world-wide.  Diffraction  tomography  affects 
the  amplitude  and  geometry  of  the  dispersion  features  primarily  at  periods  longer  than  about  100 
s.  As  documented  by  Ritzwoller  et  al.  (2002a),  diffraction  tomography,  on  average,  produces 
larger  anomalies  in  speed  in  a  period-dependent  band  of  spherical  harmonic  degrees,  and 
diffraction  and  Gaussian  tomography  maps  de-correlate  past  a  critical  spherical  harmonic  degree 
that  also  depends  on  period.  Figures  7a  and  7b  show  two  examples  of  surface  wave  dispersion 
maps  across  the  region  of  study  constructed  with  diffraction  tomography. 
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Figure  6:  Resolution  of  the  dispersion  maps  averaged  across  the  region  of  study. 


l)  Rayleigh  group  40  s  (b)  Rayleigh  group  100  s 
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Figure  7:  Selected  dispersion  maps  and  horizontal  slices  of  the  model  at  the  periods  and  depths 
indicated. 
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The  shear  velocity  model  is  constructed  using  a  Monte-Carlo  method,  which  is  described  in  detail 
by  Shapiro  &;  Ritzwoller  (2002a).  The  shear  velocity  model  of  Shapiro  &  Ritzwoller  (2002a)  was 
based  on  Gaussian  tomography  and  is  referred  to  here  as  CUB  1.0.  The  3-D  shear  velocity  model 
we  use  here  is  based  on  diffraction  tomography  and  we  will  refer  to  it  as  CUB2.0.  These  models 
only  include  vp  in  the  crust.  The  vp  part  of  CUB 2.0  in  the  mantle  will  be  derived  by  conversion 
from  vs  using  two  different  methods,  and  we  will  distinguish  between  these  two  vp  models  by 
introducing  suffices  to  the  names  (CUB2.0-TH,  CUB2.0JEMP),  as  discussed  below. 

The  inversion  for  a  velocity  profile  is  performed  at  each  node  on  a  2°  x  2°  grid  world- wide,  and 
produces  an  ensemble  of  acceptable  models  down  to  a  depth  of  400  km.  The  model  is  constrained 
by  a  variety  of  a  priori  information,  including  the  initial  crustal  model  CRUST2.0  (G.  Laske, 
personal  communication)  and  the  initial  mantle  model  S20A  of  Ekstrom  Sz  Dziewonski  (1998). 
Perturbations  are  allowed  within  specified  tolerances  to  both  vs  and  vp  in  the  crust,  to  vs  in  the 
mantle,  and  to  Moho  depth.  The  model  is  radially  anisotropic  (vsh  ^  vsv)  from  the  Moho  to  a 
variable  depth  that  averages  about  200  km.  The  strength  of  radial  anisotropy  is  constrained  to 
decrease  monotonically  with  depth  and  crustal  speeds  are  constrainted  to  increase  monotonically 
with  depth.  Shapiro  &  Ritzwoller  (2002a)  fully  describes  the  set  of  constraints  and  a  priori 
information. 

Figure  8  presents  an  example  of  the  results  of  the  Monte-Carlo  inversion  at  a  single  point.  The 
center  and  half- width  of  the  ensemble  of  acceptable  models  summarize  the  model  at  each  depth. 
The  effective  isotropic  model,  vs,  is  defined  as  the  average  of  vs^  and  vsv  at  each  depth.  Figures 
8b-d  show  that  model  uncertainties  are  reduced  and,  hence,  vertical  resolution  is  improved 
appreciably  by  simultaneously  inverting  the  group  and  phase  speed  curves.  Using  phase  velocities 
alone  produces  large  uncertainties  in  the  crust  and,  consequently,  also  in  the  upper  mantle. 
Inverting  group  velocities  alone  reduces  the  uncertainties  in  the  crust  and  uppermost  mantle  due 
to  the  measurements  at  periods  shorter  than  40  sec,  but  uncertainties  deeper  in  the  upper  mantle 
are  larger. 

When  phase  and  the  group  velocities  are  inverted  simultaneously,  however,  uncertainties  are 
significantly  reduced  at  all  depths.  In  particular,  the  simultaneous  inversion  of  broad-band  group 
and  phase  speed  data  in  the  presence  of  a  priori  constraints  on  allowable  structures  in  the  crust 
and  upper  mantle  ameliorates  the  trade-off  between  crustal  and  upper  mantle  structures  that 
plague  inversions  of  surface  waves  in  continental  areas. 
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Figure  8:  Observed  dispersion  curves  and  results  of  inversion  for  three  data  sets,  (a)  Observed 
dispersion  curves  at  a  point  in  Iran  (30°N,  60°E).  The  black  lines  are  observations  derived  from  the 
dispersion  maps  and  the  grey  shaded  zones  represent  the  range  of  dispersion  curves  from  the  the 
ensemble  of  acceptable  models  shown  in  (b).  (b)  simultaneous  inversion  of  group  and  phase  speeds 
(CUB2.0),  (c)  group  speed  alone,  and  (d)  phase  speed  alone.  The  full  width  of  the  ensemble  of 
acceptable  models  is  indicated,  where  vsh  is  dark  grey  and  vsv  is  light  grey. 
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Figure  9:  Vertical  slices  of  the  vs  model  from  CUB2.0  plotted  as  percent  deviation  from  AK135. 
Locations  of  the  profiles  are  shown  in  Figure  1. 

Figures  7c-d  and  9  present  horizontal  and  vertical  slices  of  the  vs  model  to  demonstrate  the 
nature  of  the  heterogeneities.  It  is  noteworthy  that  the  mantle  features  inferred  from  diffraction 
tomography  tend  to  have  larger  amplitudes  and  extend  deeper  than  those  from  Gaussian 
tomography.  Based  on  the  dispersion  resolution  information  given  in  Figure  6,  we  infer  that 
average  lateral  resolution  is  between  400  and  500  km  in  the  uppermost  mantle  but  degrades  with 
depth.  As  length-scales  in  the  model  approach  the  estimated  resolution,  however,  the  amplitude 
of  heterogeneity  will  tend  to  become  underestimated. 

3.2.2  Conversion  of  vs  to  vp  in  the  Mantle 

To  be  able  to  predict  Pn  and  P  travel  times,  the  3-D  5-model  CUB2.0  must  be  converted  to  a 
P-model.  There  are  two  general  approaches  to  doing  this.  The  first  is  to  use  “empirical  scaling 
relations”  that  convert  5-wave  anomalies  into  P- wave  anomalies.  The  most  successful  of  these, 
map  shear-speed  perturbations,  5vs,  relative  to  a  reference  5- model,  vso,  to  compressional-velocity 
perturbations,  Svp,  relative  to  a  reference  P-model,  vPo,  where  d  In  vp/ d  In vs  is  then  taken  to  be 
an  empirically  constrained  constant  that  may  be  a  function  of  depth,  but  is  usually  depth 
invariant.  The  second  approach  is  to  use  a  “theoretical  conversion”  based  on  laboratory 
measurements  of  thermoelastic  properties  of  mantle  minerals  and  on  models  of  the  average 
mineralogical  composition  of  the  mantle.  The  vp  model  that  results  from  CUB2.0  vs  by  the 
theoretical  conversion  will  be  referred  to  here  as  CUB2.0-TH  and  the  vp  model  derived  from  the 
empirical  scaling  relation  (d\xivs/d\nvv  =  1.9)  in  which  the  1-D  reference  model  is  AK135  will  be 
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called  CUB2.0-EMP.  We  convert  only  isotropic  vs  to  vp.  In  the  radially  anisotropic  part  of 
CUB2.0  we,  therefore,  use  (vsv  +  vsh)/2. 

We  prefer  the  theoretical  conversion  from  vs  to  vp  for  two  reasons.  First,  as  we  will  show  below, 
the  theoretical  conversion  appears  to  work  somewhat  better  in  that  the  regional  P  and  Pn 
empirical  phase  path  anomalies  are  fit  better  by  travel  times  predicted  by  CUB2.CLTH  than  by 
CUB2.0_EMP.  Second,  the  theoretical  conversion  leads  naturally  to  future  improvement.  It  can 
be  regionally  tuned  in  a  physically  meaningful  way  by  modifying  the  mineralogical  composition 
and  temperatures  within  the  anelastic  model,  and  it  can  be  updated  as  better  mineralogical  data 
become  available. 

The  theoretical  conversion  between  vs  and  vp  in  the  mantle  is  mediated  by  a  conversion  to 
temperature.  There  have  been  numerous  previous  studies  that  have  explored  the  relationship 
between  the  seismic  velocities,  temperature,  and  composition  (e.g.,  Duffy  &  Anderson,  1989; 
Graham  et  al.,  1989;  Furlong  et  ai,  1995;  Sobolev  et  a/.,  1996;  Goes  et  al.,  2000;  Rohm  et  al. , 
2000;  Trampert  et  al,  2001;  van  Wijk  et  al.,  2001).  The  method  we  use  here  is  based  on  that  of 
Goes  et  al.  (2000),  and  is  described  in  detail  by  Shapiro  Sz  Ritzwoller  (2002b).  The  mantle  is 
considered  to  be  composed  of  four  principal  minerals.  The  elastic  moduli  and  density  can  be 
computed  for  each  mineral  independently  as  a  function  of  temperature,  pressure,  and  iron 
content,  extrapolating  values  of  each  quantity  at  surface  conditions  to  depth  with  zero  iron 
content  by  using  experimentally-determined  partial  derivatives.  For  a  specific  mineralogical 
composition,  the  elastic  moduli  and  density,  and  hence  the  seismic  velocities,  are  computed  using 
the  Voigt-Reuss-Hill  mixing  scheme.  The  result  holds  at  high  frequencies  or  very  low 
temperatures  where  anelasticity  contributes  minimally.  At  mantle  temperatures  and  seismic 
frequencies,  however,  the  temperature-velocity  relation  must  include  a  correction  for  physical 
dispersion.  Because  upper-most  mantle  Q  is  poorly  known,  we  use  the  attenuation  model  of 
Minster  k  Anderson  (1981)  that  relates  Q  to  temperature.  We  follow  Sobolev  et  al.  (1996)  in 
specifying  the  parameters  in  this  conversion  (frequency  exponent,  activation  energy  and  volume), 
but  we  calibrate  the  amplitude  of  the  correction  based  on  the  average  shear-velocity  in  the  region 
of  study  and  an  assumed  average  upper  mantle  temperature  of  1400°  C  at  200  km  depth.  We  use 
a  single  mineralogical  composition  here  for  the  entire  region  of  study,  the  average  off-cratonic 
continental  composition  advocated  by  McDonough  k  Rudnick  (1998):  68%  Olivine,  18% 
Orthopyroxene,  11%  Clinopyroxene,  and  3%  Garnet  with  an  Iron: Magnesium  ratio  of  10%. 
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Figure  10:  vp  models  at  the  midpoint  of  two  profiles,  (a)  Grid  used  for  the  theoretical  vs  to  vp 
conversion,  (b)  The  theoretical  conversion  presented  as  the  logarithmic  scaling  factor  dlnvs/dlnvp. 
The  solid  black  line  is  vs  from  the  1-D  model  AK135  and  the  dashed  line  is  vs  converted  from  AK135 
vp  by  the  theoretical  conversion. 

Figure  10a  shows  the  resulting  vs  to  vp  theoretical  conversion.  Figure  10b  displays  this  conversion 
presented  as  the  logarithmic  scaling  relation,  dlnvs/dlnvp:  which  varies  with  both  vs  and  depth. 
The  vs  profile  from  AK135  is  overplotted,  nearly  paralleling  the  contours  of  the  theoretical 
predictions.  This  illustrates  why  depth-independent  values  of  the  scaling  relation  tend  to  work 
fairly  well  in  the  upper  mantle.  For  the  values  of  vs  in  AK135,  the  theoretical  prediction  for  the 
scaling  relation  is  dlnvs/d\nvp  ~  1.6  —  1.8.  Figure  10b  also  shows  that  the  vs  profile  converted 
from  the  AK135  vv  profile  by  the  theoretical  conversion  agrees  fairly  well  with  the  vs  profile  in 
AK135  at  depths  below  about  100  km.  The  theoretical  conversion  between  vs  and  vp  differs 
appreciably  from  the  vs  and  vp  parts  of  AK135  above  about  100  km. 

Thus,  in  the  shallower  parts  of  the  mantle,  vp  computed  using  the  theoretical  conversion  will  differ 
appreciably  from  vp  computed  using  the  logarithmic  scaling  relation  applied  to  reference  values 
from  the  1-D  model  AK135.  Figure  11  exemplifies  this  by  showing  representative  vp  profiles  from 
tectonic  (e.g.,  A  —  A1)  and  platform  (e.g.,  B  —  B! )  regions  of  CUB2.0-TH  and  CUB2.0JEMP. 
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Figure  11:  vp  model  at  the  mid-points  of  profiles  (a)  A  —  A!  and  (b)  B  —  Bf  from  Figures  1  and  9. 
The  solid  line  represents  the  theoretical  conversion  from  vs  (i.e.,  CUB2.CLTH)  and  the  dashed  line 
is  from  the  empirical  scaling  relation  dlnvs/d\nvp  =  2.0  in  which  AK135  is  used  as  the  reference 
(i.e.,  CUB2.0JEMP). 

3.2.3  Travel  Times  and  Correction  Surfaces 

To  compute  travel  times  for  Pn  and  P  we  use  a  modified  version  of  the  2-D  ray  tracer  developed 
by  Cerveny  &  Psencik  (1984).  Villasenor  et  al.  (2002)  shows  that  this  code  applied  to  the  CUB 
3-D  models  agrees  well  with  travel  times  from  the  finite  difference  code  of  Podvin  &  Lecompte 
(1991)  modified  to  be  applied  on  a  spherical  earth.  They  also  show  that  travel  time  variations 
from  3-D  excursions  in  regional  rays  through  the  CUB  3-D  models  are  small  and  do  not 
substantiate  the  computational  cost  of  computing  3-D  rays.  Even  in  2-D,  Pn  is  a  very  complex 
phase  as  the  rays  and  travel  time  curves  shown  in  Figure  12  illustrate. 

The  travel  times  of  the  Pn  phase  from  the  3-D  P-models  derived  with  the  two  conversion  schemes 
can  differ  appreciably,  as  Figure  13  shows.  In  a  geographically  averaged  sense,  CUB2.CLEMP 
produces  Pn  travel  times  that  are  about  1  sec  faster  than  Pn  times  from  CUB2.0-TH.  This 
discrepancy  is  strongest  in  tectonic  areas,  but  in  shield  areas  the  theoretical  conversion  actually 
tends  to  produce  faster  times  at  distance  beyond  about  100  km  because  the  deeper  parts  of  the 
model  beneath  the  shallow  mantle  lid  are  faster  and  rays  tend  to  dive  deeper.  The  difference 
between  the  empirical  and  theoretically  converted  Pn  travel  times  is  largest  at  distances  beyond 
about  2000  km.  This  is  probably  because  second-order  pressure  dependencies  are  needed  in  the 
theoretical  conversion  for  rays  turning  deeper  than  about  250  km.  Future  advances  in  the 
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Figure  12:  Rays  traced  through  the  vertical  slices  of  the  3-D  model  CUB2.CLTH  and  corresponding 
travel  time  curves,  (a)  and  (b)  Pn  rays  traced  through  the  vertical  slices  of  CUB2.CLTH  from 
Figures  1  and  9.  P-velocities  in  the  upper  mantle  underlie  the  rays,  (c)  and  (d)  The  reduced  Pn 
travel  time  curves  for  each  slice  (relative  to  8.0  km/s),  shown  as  solid  lines,  are  compared  with  the 
Pn  first  arrivals  for  AK135.  This  demonstrates  the  complexity  of  Pn. 


Pn  time  difference  (s) 


Figure  13:  Effect  of  the  vs  to  vp  conversion  on  Pn  travel  times,  (a)  Difference  between  Pn  travel 
times  from  CUB2.CLTH  and  CUB2.0-EMP  centered  on  stations  KIV  (Kislovodsk,  Russia).  Positive 
values  mean  CUB2.0_TH  travel  times  are  faster  than  those  from  CUB2.CLEMP.  (b)  Solid  line  shows 
the  Pn  travel  time  difference  for  profile  B  —  Br  and  the  dashed  line  for  A  —  A!  from  Figures  1  and 
9.  One  standard  deviation  errors  bars  show  the  geographical  average  across  the  region  of  study 
(Figure  1). 
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Figure  14:  Example  of  Pn  and  P  travel  time  correction  surfaces  computed  for  a  surface  source  for 
station  KIV.  The  corrections  here  are  relative  to  the  1-D  model  IASP91. 

theoretical  conversion  based  on  finite-strain  theory  may  correct  this  effect,  but  for  the  present 
study  it  is  a  moot  point,  at  least  for  Pn,  which  is  rarely  observed  beyond  2000  km.  The  use  of  the 
theoretical  conversion  for  longer  baseline  P  phases,  however,  will  need  to  take  this  into 
consideration. 

For  the  validation  exercises  in  later  sections,  we  compute  “travel  time  correction  surfaces”  for  Pn 
and  P  for  all  of  the  854  stations  in  the  groomed  arrival  time  data  set.  Separate  surfaces  are 
computed  for  a  discrete  set  of  hypothesized  event  depths  (0,  10,  20,  30,  40  km)  and  are  then 
interpolated  to  provide  correction  surfaces  for  events  between  this  set  of  depths.  These  travel  time 
surfaces  depict  the  predicted  travel  times  from  events  on  a  grid  of  epicentral  locations  and  depths 
observed  at  a  particular  station,  and  are  presented  relative  to  the  travel  time  predicted  from  a  1-D 
model.  Figure  14  presents  examples  of  regional  station-centered  travel  time  correction  surfaces  for 
P  and  Pn  computed  from  CUB2.0_TH  relative  to  the  travel  times  from  IASP91.  Figure  15 
presents  examples  of  event- centered  correction  surfaces  relative  to  travel  times  from  AK135. 
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travel  time  correction  relative  to  AK1 35  (s) 


Figure  15:  Cluster-centered  travel  time  correction  surfaces  for  four  cluster  regions.  The  colored 
contours  are  the  predictions  from  the  3-D  model  CUB2.CLTH  referenced  to  the  travel  time  from 
the  1-D  model  AK135.  The  symbols  are  the  empirical  phase  path  anomalies,  color-coded  similarly 
to  the  model  predictions.  The  cluster  time  baseline  shift  is  indicated  in  Tables  5  and  6.  (a)  and  (b) 
are  explosions  in  the  Azgir  and  Lop  Nor  clusters  and  (c)  and  (d)  are  earthquakes  in  the  Chamoli 
and  Racha  clusters,  respectively. 
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3.3  Fit  to  Regional  Travel  Times 


The  first  assessment  of  the  3-D  model  is  to  determine  how  well  regional  travel  times  predicted  by 
the  3-D  P-wave  models  CUB2.0-TH  and  CUB2.CLEMP  fit  well  determined  regional  travel  time 
observations.  This  test  is  made  possible  by  the  empirical  phase  path  anomaly  data  set  described 
in  section  2.  We  compute  travel  times  from  each  model  tested  by  using  the  GT  locations,  depths, 
and  origin  times.  We  will  argue  that  the  systematics  of  misfit  establish  that  the  3-D  model 
greatly  improves  the  fit  to  regional  travel  times  over  1-D  models  and  that  the  theoretical 
conversion  from  vs  is  preferable  to  the  empirical  scaling  relation. 

We  perform  comparisons  both  cluster-by-cluster  and  aggregated  over  clusters,  segregating  the 
results  by  event  type  (e.g.,  explosions,  GTS  or  GT10  earthquakes).  Examples  of  cluster-specific 
comparisons  are  presented  in  Figures  15  and  16.  In  each  case,  we  allow  an  offset  to  the  predicted 
correction  surface  to  minimize  the  rms  travel  time  residual.  A  similar  shift  was  also  introduced  in 
the  HDC  analysis  and  in  the  construction  of  the  empirical  phase  path  anomalies  discussed  in 
section  2,  referred  to  as  “cluster  time  baseline  shifts”. 

These  shifts  are  listed  in  Table  5  in  which  the  1-D  model  AK135  and  both  regional  and 
teleseismic  data  were  used.  Using  just  regional  data,  they  are  listed  in  Tables  5  -  7  for  the  3-D 
model  CUB2.CLTH  and  the  1-D  model  AK135.  Tables  5  -  7  also  provide  information  about  misfit 
and  correlation  between  the  model  predicted  travel  times  and  the  empirical  phase  path  anomalies. 
Note  that  correlation  is  defined  with  respect  to  perturbations  relative  to  AK135.  The  correlation 
of  AK135  travel  times  with  the  empirical  phase  path  anomalies  is  zero,  therefore,  because  the 
perturbations  are  all  zero.  For  this  reason,  the  correlations  for  AK135  are  left  out  of  Tables  5  -  7. 
Figure  15  presents  examples  of  cluster-centered  travel  time  correction  surfaces  produced  from 
CU2.CLTH  with  the  empirical  path  anomalies  for  Pn  overlain.  Plots  such  as  these  establish  the 
qualitative  agreement  between  the  empirical  path  anomalies  and  the  model  predicted  travel  times, 
at  least  with  respect  to  the  large-scale  features  of  the  correction  surfaces.  More  quantitative 
comparisons  are  seen  in  Figure  16,  which  also  demonstrates  the  effect  of  phase  re-identification 
using  the  3-D  model.  Aggregated  comparisons,  such  those  shown  in  Figure  17  and  the  averages 
presented  in  Tables  5-7,  demonstrate  that  the  model  fit  to  the  empirical  path  anomalies  depends 
strongly  on  the  GT  level,  with  explosions  being  fit  best,  followed  by  GT5  and  finally  GT10 
earthquakes.  The  rms-misfit  averaged  over  the  explosion  clusters  for  CUB2.CLTH  is  about  1.04 
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Figure  16:  Comparison  between  the  P  and  Pn  travel  times  predicted  by  the  3-D  model  CUB2.0-TH 
and  the  empirical  phase  path  anomalies  for  the  four  event  clusters  shown  in  Figure  15.  The  empirical 
path  anomalies  are  on  the  horizontal  axis,  paired  with  the  model  predictions  on  the  vertical  axis. 
Circles  are  Pn  and  triangles  are  P.  Phases  have  been  re-identified  using  the  3-D  model;  bars  link 
phases  originally  identified  with  the  1-D  model  AK135  (shaded  gray)  to  the  re-identified  phases 
(unshaded).  The  cluster  time  baseline  shifts  listed  in  Tables  5  and  6  have  been  applied.  The  dashed 
lines  mark  the  location  of  ±3  sec  residuals. 
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sec,  and  misfit  averaged  over  the  GT5  and  GT10  earthquake  clusters  averages  1.24  sec  and  1.54 
sec,  respectively. 


Table  5.  Fit  to  empirica 


phase  path  anoma 


CUB2.CLTH 


Cluster  No. 

Name 

RMSa 

Correl6 

At0c 

RMSa 

A£0c 

1 

Azgir 

0.96 

0.87 

1.17 

2.63 

2.52 

2 

Balapan 

0.99 

0.74 

1.69 

1.33 

2.06 

3 

Degelen 

0.99 

0.67 

0.54 

1.09 

1.47 

4 

Lop  Nor 

1.50 

0.62 

1.66 

1.85 

1.25 

5 

Novaya 

0.73 

0.82 

0.15 

1.59 

3.14 

average^ 

:  1.04 

0.75 

1.04 

1.70 

1  2.09 

ies  for  explosions. 


AK135 


a  Rms-misfit  between  observed  and  predicted  travel  times,  in  sec. 
b  Correlation  between  observed  and  predicted  travel  times. 
c  Cluster  time  baseline  shift  applied  to  regional  data  alone,  in  sec. 
d  Average  of  absolute  values  of  cluster  values. 


Table  6.  Fit  to  empirical  phase  path  anomalies  for  GT5  earthquakes. 


Terminology  similar  to  Table  5. 


CUB2.0-TH 

AK135 

Cluster  No. 

Name 

RMS 

Correl 

o 

-t-J 

<1 

RMS 

At0 

7 

Adana 

1.07 

0.50 

-1.76 

1.18 

-2.74 

8 

Bhuj 

1.37 

0.76 

-0.78 

2.00 

0.98 

9 

Chamoli 

2.10 

0.65 

-0.49 

2.36 

-0.90 

10 

Duzce 

1.32 

0.52 

-0.14 

1.38 

-0.90 

11 

Erzin 

1.22 

0.77 

-1.78 

2.14 

-2.91 

12 

Garm 

1.17 

0.86 

-0.89 

2.56 

-1.71 

13 

Hoceima 

0.76 

0.43 

-0.67 

0.80 

-0.33 

14 

Izmit 

1.04 

0.57 

-0.10 

1.14 

-0.83 

15 

Koyna 

1.02 

0.85 

-0.94 

1.62 

2.31 

16 

Racha 

1.39 

0.71 

-0.15 

1.89 

-1.25 

17 

Siberia 

1.13 

0.69 

-0.65 

1.55 

-0.98 

average 

1.24 

0.67 

0.76 

1.70 

1.44 
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Figure  17:  Overall  comparison  between  the  P  and  Pn  travel  times  predicted  by  the  3-D  model 
CUB2.0-TH  and  the  empirical  phase  path  anomalies.  The  cluster  time  baseline  shifts  listed  in 
Tables  5-7.  have  been  applied  and  misfit  statistics  are  also  listed  in  these  tables. 
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Table  7.  Fit  to  empirical  phase  path  anomalies  for  GT10  earthquakes. 


Terminology  similar  to  Table  5. 


CUB2.0-TH 

AK135 

Cluster  No. 

Name 

RMS 

Correl 

A  t0 

RMS 

At0 

19 

Aqaba 

1.60 

0.31 

-0.11 

1.54 

-0.74 

20 

Gubal 

1.02 

0.76 

0.09 

1.78 

-0.98 

21 

Jiashi 

2.07 

0.77 

0.40 

2.99 

-1.57 

22 

Spitak 

1.24 

0.84 

-0.08 

2.28 

-1.45 

23 

Tab  as 

1.74 

0.37 

-0.71 

1.88 

-2.42 

average 

1.54 

0.61 

0.28 

2.10 

1.43 

These  values  are  to  be  contrasted  with  misfit  using  the  1-D  model  AK135:  1.70  sec,  1.70  sec,  and 
2.10  sec,  respectively.  The  3-D  model,  therefore,  delivers  a  64%  variance  reduction  relative  to  the 
1-D  model  for  the  explosions  and  a  46%  variance  reduction  for  the  earthquakes.  The  fact  that  the 
GT  level  affects  the  fit  to  the  regional  travel  times  is  not  surprising,  but  does  lend  support  to  the 
procedure  that  established  the  GT-level  discussed  in  section  2.  Figures  16  and  17  display  as 
dashed  lines  the  ±3  sec  residual  levels  to  demonstrate  that  misfits  at  this  level  are  rare,  especially 
for  explosions. 

In  the  location  exercise  below,  we  will,  therefore,  use  ±3  sec  as  the  cut-off  for  out-lier 
identification.  Residuals  with  absolute  values  greater  than  3  sec  are  suspicious  measurements  or 
are  phases  that  are  particularly  difficult  to  identify.  Tables  5  -  7  also  list  the  cluster  time  baseline 
shifts  for  CUB2.CLTH  and  AK135  using  the  regional  empirical  phase  path  anomalies  alone  (i.e.,  in 
contrast  to  Table  3  in  which  teleseismic  travel  times  were  used).  The  cluster  time  baseline  shifts 
are  introduced  to  compensate  for  origin  time  errors  in  the  GT  data  base  and  to  aid  in  phase 
identification,  but  they  also  act  to  correct  model  errors  in  both  the  crust  and  uppermost  mantle, 
particularly  beneath  the  cluster.  Tables  5-7  show  that,  on  average,  the  cluster  time  baseline 
shifts  are  systematically  smaller  for  CUB2.CLTH  than  for  AK135  (explosions:  1.04  sec  versus  2.09 
sec;  GT5  earthquakes:  0.76  sec  versus  1.44  sec;  GT10  earthquakes:  0.28  sec  versus  1.43  sec).  In 
fact,  they  are  smaller  for  all  but  two  of  the  clusters  which  have  anomalously  small  shifts  for  the 
1-D  model.  This  is  presumably  because  the  3-D  model  CUB2.0_TH  more  accurately  models 
near-cluster  structure  than  the  1-D  model  AK135.  (Results  are  not  presented  for  the  Aden  or 
Sahara  clusters  because  there  are  too  few  regional  data  for  meaningful  statistics.) 
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Because  the  cluster  time  baseline  shifts  may  obscure  some  problems  with  the  model,  it  is  also 
important  to  perform  travel  time  comparisons  without  the  shifts.  Figure  18  compares  the 
empirical  P  and  Pn  path  anomalies  with  the  travel  times  predicted  by  a  number  number  of 
models  (1-D:  AK135,  PREM;  3-D:  CU2.CLTH,  CU2.0.EMP,  AK135+CRUST5.1+S20A_TH)  as  a 
function  of  epicentral  distance.  The  model  AK135+CRUST5.1+S20A_TH  is  a  vp  model  of  the 
crust  and  mantle  constructed  by  placing  CRUST5.1  (Mooney  et  ai,  1998)  on  top  of  P- wave 
speeds  converted  from  the  3-D  mantle  vs  model  S20A  (Ekstrom  &  Dziewonski,  1998).  The  model 
S20A  was,  in  fact,  developed  as  an  aspherical  perturbation  to  PREM,  but  because  PREM  was  not 
developed  to  fit  regional  travel  times,  we  replace  it  with  AK135  as  the  spherical  reference  for  this 
model.  Finally,  we  use  the  theoretical  conversion  from  vs  to  vp  to  compute  the  vp  part  of  the 
mantle  model;  hence  the  suffix  on  the  name  AK135+CRUST5.1+S20A_TH. 

Figures  18a  -  18e  show  residuals  for  every  empirical  path  anomaly  in  the  data  base.  These 
residuals  are  averaged  in  distance  bins  and  broken  into  a  vertical  offset  (or  bias)  and  the  standard 
deviation  around  the  offset  as  shown  in  Figures  18f  and  18g.  We  divide  the  comparison  in  this 
way  because  it  is  possible  for  a  model  to  fit  the  geographical  pattern  of  travel  time  residuals  well 
but  be  systematically  biased.  The  overall  rms-residual  would,  therefore,  be  large,  but  that  alone 
would  provide  little  insight  into  the  reason  behind  the  large  residuals. 

Figures  18f  and  18g  demonstrate  that  the  1-D  model  AK135,  which  was  constructed  to  fit  regional 
travel  times  on  average,  fits  the  empirical  path  anomalies  better  than  PREM,  which  displays  a 
bias  of  about  4  sec  at  distances  less  than  1500  km.  This  is  largely  due  to  thin  crust  in  PREM 
(~20  km  after  the  ocean  was  removed  for  continental  application).  The  travel  times  from  AK135 
are  also  biased,  but  at  a  lower  level  than  PREM,  and  the  distance  trend  is  less  severe.  The 
pattern  of  bias  for  CUB2.0-EMP  reflects  the  bias  in  AK135,  probably  because  of  the  reliance  of 
the  empirical  vs  to  vp  conversion  on  AK135.  The  travel  times  from  the  two  3-D  models  that  have 
been  converted  to  vp  with  the  theoretical  conversion  (CUB2.CLTH,  AK135+CRUST5.1+S20A) 
are  much  less  biased.  This  is  one  of  the  reasons  for  our  preference  for  the  theoretical  vs  to  vp 
conversion.  We  will,  therefore,  not  consider  CUB2.0JEMP  in  the  location  exercises  below. 

The  fit  to  the  geographical  pattern  of  empirical  phase  path  anomalies  is  revealed  in  the  standard 
deviation  shown  in  Figure  18g.  The  3-D  models  all  fit  better  than  the  1-D  models.  CUB2.CLTH 
has  only  a  slightly  smaller  standard  deviation  than  CUB2.CLEMP.  CUB2.CLTH  fits  the  pattern  of 
empirical  phase  path  anomalies  somewhat  better  than  AK135+CRUST5.1+S20A.  This  is  a 
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Figure  18:  Residuals  between  the  empirical  phase  path  anomalies  for  Pn  and  P  and  predictions. 
(a)-(e)  Residuals  from  a  number  of  models,  plotted  versus  epicentral  distance.  Running  average 
(bias)  and  standard  deviation  are  shown  as  the  black  lines.  Cluster  time  baseline  shifts  have  not 
been  applied,  (f)  -  (g)  The  running  residual  bias  and  standard  deviation  shown  for  the  various 
models. 
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promising  because  CUB2.CLTH  was  built  by  perturbing  AK135+CRUST5.1+S20A.  This 
perturbation  was  accomplished  by  introducing  many  regional  group  speed  measurements  in  the 
attempt  to  improve  lateral  and,  more  importantly,  vertical  resolution  in  the  resulting  model  (e.g., 
Figure  8).  There  was  no  guarantee  that  this  procedure  would  produce  a  vp  model  that  would 
improve  the  fit  to  regional  P-phases,  but  these  results  establish  that  it  has.  This  points  the  way  to 
future  advancements  in  vs  models  providing  further  improvements  in  regional  location  capabilities. 

3.4  Regional  Location  Experiment 

The  location  experiment  is  performed  using  only  P  and  Pn  travel  times  from  the  groomed  arrival 
time  data  set,  observed  within  20°  of  the  epicenter.  The  figure-of-merit  is  the  ability  to  locate 
those  events  in  the  GT  data  set  with  location  accuracy  of  GT5  or  better.  We  compare  the 
location  capabilities  of  the  3-D  model  CUB2.0-TH  to  that  of  the  1-D  model  AK135.  Examples  of 
locations  using  these  models  are  shown  in  Figure  19  for  the  Lop  Nor  and  Racha  clusters. 
Qualitatively,  it  can  be  seen  that  the  3-D  model  improves  the  agreement  with  the  GT  locations 
and  reduces  systematic  bias. 

The  geometry  of  the  network  of  reporting  stations  plays  a  major  role  in  location  accuracy.  Figure 
20  presents  examples  of  how  location  error  varies  with  the  number  of  stations  and  open  azimuth. 
These  results  are  determined  by  performing  locations  repeatedly  using  randomly  chosen  subsets  of 
the  reported  data  in  the  groomed  arrival  time  data  set  with  a  specified  number  of  recording 
stations. 

The  location  systematics  for  the  3-D  model  are  clear  and  reasonable:  as  the  number  of  stations 
decreases  and  open  azimuth  increases  the  location  accuracy  degrades.  We  aim,  however,  for  the 
statistics  on  location  accuracy  to  reflect  the  capabilities  of  the  model  rather  than  the  vagaries  of 
station  geometry.  It  is,  therefore,  important  to  limit  the  role  that  variations  in  network  geometry 
between  different  event  regions  play  on  the  reported  locations.  For  this  reason,  the  location 
statistics  reported  in  Tables  8-10  are  only  for  the  GT  events  with  open  azimuths  less  than  180° 
for  regional  P  and  Pn  data.  This  reduces  the  GT  data  set  by  more  than  60%,  from  989  to  366 
events.  The  number  of  explosions  is  most  severely  reduced,  from  312  to  38  events,  but  the  number 
of  events  remains  high  enough  to  draw  statistically  meaningful  conclusions.  Using  the  3-D  model, 
the  location  accuracy  improves  systematically  with  the  confidence  in  the  GT  location:  5.1  km,  7.2 
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★  GT  location  a  3D  location  •  ID  location 


Figure  19:  Mislocation  vectors  for  two  event  clusters  for  (left  column)  explosions  on  the  Lop  Nor 
test  site  and  (right  column)  earthquakes  near  Racha.  Georgia,  (a)  and  (b)  3-D  model  (CUB2.CLTH) 
locations,  (c)  and  (d)  1-D  model  (AK135)  locations.  Stars  mark  GT  locations,  triangles  and  circles 
are  3-D  and  1-D  model  locations,  respectively. 
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Figure  20:  Mislocation  using  the  3-D  model  CUB2.0-TH,  presented  as  a  function  of  number  of 
stations  and  open  azimuth,  determined  by  randomly  choosing  subsets  of  the  reported  travel  times 
for  four  events:  (a)  Lop  Nor,  May  26,  1990;  (b)  Lop  Nor,  Sept  25,  1992;  (c)  Racha,  May  3,  1991; 
(d)  Racha,  July  4,  1991. 


km,  and  12.3  km  errors  for  explosions  (GT0-GT2),  GT5  earthquakes,  and  GT10  earthquakes, 
respectively.  By  the  following  reasoning,  we  believe  that  these  results  are  consistent  with  an 
intrinsic  location  accuracy  of  about  5  km  for  the  3-D  model. 


Table  8.  Location  results  for  nuclear  explosions. 


CUB2.0-TH 

AK135 

Cluster  No. 

Name 

Events0 

Stations5 

Error6 

! 

RMSd 

Error6  : 

RMSd 

1 

Azgir 

6 

19 

4.2 

1.08 

24.0 

1.49 

2 

Balapan 

10 

25 

3.9 

1.10 

11.2 

1.22 

3 

Degelen 

3 

18 

3.1 

1.01 

12.0 

1.36 

4 

Lop  Nor 

16 

45 

6.1 

1.22 

9.6 

1.49 

5 

Novaya 

3 

15 

!  8.0 

1.03 

13.6 

1.57 

overall6 

38 

24 

5.1 

1.09 

14.1 

1.43 

a  Number  of  events  with  open  azimuth  <  180°  and  >  5  stations. 
b  Median  number  of  stations. 

c  Average  distance  from  estimated  to  GT  location,  in  km. 
d  Rms  difference  between  observed  and  predicted  travel  times,  in  sec. 
e  Total  or  average  of  cluster  values. 
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Table  9.  Location  results  for  GTS  earthquakes. 


Terminology  similar  to  Table  8. 


CUB2.0.TH 

AK135 

Cluster  No. 

Name 

Events 

Stations 

Error 

RMS 

Error 

RMS 

7 

Adana 

19 

23 

00 

bo 

1.04 

8.1 

1.17 

8 

Bhuj 

4 

15 

11.5 

1.00 

8.5 

1.42 

9 

Chamoli 

50 

17 

10.8  . 

1.20 

20.8 

1.54 

10 

Duzce 

21 

45 

7.3 

1.23 

6.6 

i 

1.33 

11 

Erzin 

8 

33 

6.0 

1.29 

9.5 

1.61 

12 

Garm 

26 

21 

5.2 

1.12 

9.7 

2.17 

13 

Hoceima 

21 

26 

4.4 

0.90 

7.5 

0.86 

14 

Izmit 

8 

116 

5.1 

1.07 

5.3 

1.10 

15 

Koyna 

8 

10 

5.5 

1.13 

10.7 

1.23 

16 

Racha 

34 

20 

8.2 

1.15 

21.8 

1.43 

17 

Siberia 

8 

33 

6.3 

!  1.13 

11.0 

1.31 

overall 

207 

32 

7.2 

1.12 

10.9 

1.38 

Table  10.  Location  results  for  GT10  earthquakes. 


Terminology  similar  to  Table  8. 


CUB2.0-TH 

AK135 

Cluster  No. 

Name 

Events 

Stations 

Error  j 

RMS 

Error 

RMS 

19 

Aqaba 

14 

47 

17.6 

1.17 

13.8 

1.27 

20 

Gubal 

10 

16 

9.8 

1.12 

14.9 

1.31 

21 

Jiashi 

61 

12 

16.6 

0.87 

20.0 

1.69 

22 

Spitak 

10 

32 

8.3 

1.24 

9.1 

1.69 

23 

Tab  as 

26 

19 

9.3 

1.17 

8.7 

1.27 

overall 

121 

25 

12.3 

1.11 

13.3 

1.44 

Assuming  that  the  regional  location  errors  and  the  reported  GT  confidence  levels  are 
uncorrelated,  we  would  expect  a  location  error  of  about  (52  +  52)1/2  km  ~  7  km  for  GTS  events  if 
the  intrinsic  location  error  is  itself  5  km,  and  (102  +  52)1/2  km  ~  12  km  for  GT10  events.  These 
expectations  are  very  similar  to  to  the  results  of  our  tests,  and  give  us  further  confidence  in  the 
GT  levels  reported  in  Table  3. 

Location  errors  for  the  1-D  model  do  not  trend  as  simply  with  GT  level.  Location  errors  are  14.1 
km,  10.9  km,  and  13.3  km,  respectively,  for  explosions,  GT5  earthquakes,  and  GT10  earthquakes. 
The  location  error  for  the  explosions  is  elevated  by  the  difficulties  that  the  1-D  model  has  in 
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locating  the  Azgir  events.  The  reason  is  the  nature  of  the  3-D  structure  near  Azgir:  fast 
continental  platforms  to  the  north  and  slow  tectonic  regions  to  the  south.  The  3-D  model  corrects 
for  this  large-scale  variability  very  well,  but  the  1-D  model  produces  systematically  biased 
locations  similar  to  its  performance  at  Racha  (Figure  19d).  It  is,  nevertheless,  reasonable  to 
conclude  that  the  1-D  model  possesses  an  average  intrinsic  location  accuracy  of  about  10  km  from 
which  we  would  expect  estimated  location  errors  of  about  10  km,  11  km,  and  14  km,  respectively, 
for  explosions,  GT5  earthquakes,  and  GT10  earthquakes.  This  accuracy  is  much  more 
geographically  variable  than  for  the  3-D  model,  however. 

These  estimates  of  intrinsic  location  accuracies  (5  km  for  the  3-D  model,  10  km  for  the  1-D 
model)  have  been  determined  for  relatively  large  events  for  which  a  large  number  of  regional 
phases  were  available.  It  is  important  in  a  number  of  applications,  particularly  in  nuclear 
monitoring,  to  understand  how  location  accuracy  degrades  as  the  network  of  recording  stations 
becomes  increasingly  sparse.  Figure  21  addresses  this  question  by  considering  location  accuracy 
as  a  function  of  random  subsets  of  the  reported  stations.  The  approach  is  similar  to  that  taken  in 
Figure  20,  but  here  we  have  aggregated  the  results  over  different  types  of  events:  explosions  in 
Figure  21a  and  GTS  events  in  Figure  21b. 

The  averages  and  standard  deviations  of  the  distributions  of  locations  are  presented  as  a  function 
of  the  number  of  stations  retained,  with  the  constraint  that  open  azimuth  remains  less  than  180°. 
Note  that  explosion  locations  for  the  3-D  model  degrade  from  about  5  km  using  a  large  number  of 
stations  to  about  8  km  on  average  for  5  reported  stations.  GTS  earthquake  locations  using  the 
3-D  model  similarly  degrade  from  averages  of  about  11  km  to  17  km.  Figures  21a  and  21b  also 
show  that  (1)  the  average  location  error  grows  as  the  number  of  stations  decreases  with 
degradation  setting  on  at  about  10  stations,  (2)  the  mean  and  standard  deviation  of  the 
distributions  are  smaller  for  the  3-D  model  than  for  the  1-D  model,  and  (3)  the  mean  and 
standard  deviation  are  smaller  for  explosions  than  for  GTS  earthquakes. 

What  is  not  shown  in  Figures  21a  and  21b  is  that  the  distributions  of  location  error  for  the  two 
models  are  highly  correlated.  Therefore,  even  though  the  location  error  distributions  for  the  two 
models  for  the  GTS  earthquakes  overlap  appreciably,  the  3-D  model  location  is  better  than  the 
1-D  model  location  in  no  less  than  70%  of  the  cases  considered  as  21c  shows.  For  the  GTS 
earthquakes,  this  number  remains  approximately  constant  with  the  number  of  reporting  stations. 
For  the  explosions,  however,  it  degrades  from  more  than  90%  for  large  numbers  of  stations  used 
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Figure  21:  Location  error  for  the  3-D  and  1-D  models,  plotted  as  a  function  of  the  number  of 
reporting  stations,  determined  by  randomly  choosing  subsets  of  the  reported  travel  times  in  which 
open-azimuth  is  constrained  to  be  less  than  180°.  (a)  and  (b)  Average  (solid  lines)  and  standard 
deviation  (dashed  lines)  are  shown  for  explosions  and  GTS  earthquakes  for  CUB2.0-TH  (black 
lines)  and  AK135  (grey  lines),  (c)  Percent  of  the  events  for  which  the  location  error  is  smaller 
using  the  3-D  model  than  the  1-D  model  for  explosions  (solid  line)  and  GTS  earthquakes  (dashed 
line) , 
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to  about  75  -  80%  for  sparse  networks  containing  less  than  10  regional  stations. 


4.  Conclusions 

The  validation  tests  reported  above  yield  a  number  of  general  lessons  that  we  discuss  further  here. 

4.1  GT  Data  Bases 

Data  bases  of  earthquake  and  explosion  locations  with  high  accuracy  and  quantified  uncertainties 
(Ground  Truth  data  bases)  provide  valuable  information  needed  to  test  earth  models  and  to 
determine  seismic  location  capabilities.  Such  data  bases  are  very  difficult  to  assemble,  however,  as 
they  require  substantial  efforts  to  validate.  The  Hypocentroidal  Decomposition  (HCD)  method  of 
multiple  event  location,  in  particular,  is  very  well  suited  to  the  requirements  of  Ground  Truth 
(GT)  validation  exercises.  The  empirical  phase  path  anomalies  that  can  be  constructed  from 
clusters  of  events  provide  important  ancillary  information  about  the  capabilities  of  earth  models. 

It  is  a  popular  conception  that  compilations  of  regional  empirical  phase  path  anomalies  may  be 
used  directly  to  construct  empirical  “correction  surfaces’5  to  be  used  in  location  procedures.  Our 
experience  is  that  information  of  GT10  quality  or  better  results  only  after  strenuous  efforts,  and 
accrues  only  in  isolated  regions  of  the  globe  where  reference  events  are  available.  It  is  not  likely, 
therefore,  that  empirically-derived  regional  correction  surfaces  can  be  developed  from  GT10  data 
bases  in  a  general  sense.  It  remains  to  be  determined  if  empirical  correction  surfaces  that  derive 
from  lower  quality  locations  wTill  be  able  to  perform  as  well  as  existing  3-D  earth  models,  such  as 
CUB2.CLTH,  over  large,  diverse,  and  widely  dispersed  regions.  It  is,  for  example,  still  not  fully 
understood  how  best  to  merge  model  predictions  with  empirical  path  corrections.  We  have  shown 
here  that  high  quality  baseline  information  from  existing  3-D  models  exists  already,  and  provides 
a  context  into  which  empirical  information  of  similar  quality  can  be  introduced  in  the  future. 

4.2  S- Wave  Models 

The  3-D  models  tested  here  (CUB2.0,  S20A)  are,  in  essence,  vs  models  of  the  upper  mantle  that- 
have  been  derived  largely  from  information  about  surface  wave  dispersion.  Surface  waves  possess 
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the  salutary  characteristic  that  they  sample  the  entire  earth,  and  over  large  regions  data  coverage 
is  dense  and  relatively  homogeneous.  Another  advantage  is  that  broad-band  surface  wave 
dispersion  constrains  the  vertical  velocity  gradient  in  the  upper  mantle.  The  vertical  gradient 
controls  ray  turning  depths  which,  in  turn,  largely  control  travel  times. 

It  is  much  more  difficult  to  constrain  3-D  models  of  vp  in  the  uppermost  mantle  directly,  except  in 
rare  regions  with  exceptional  station  coverage  and  seismicity.  Pn  tomography  produces  2-D 
surfaces  which  provide  no  information  about  the  vertical  velocity  gradient  generally  needed  to 
predict  regional  travel  times.  Pn  tomographic  maps  have  been  shown  by  Ritzwoller  et  al  (2002b) 
to  correlate  well  with  uppermost  mantle  P-wave  velocities  computed  from  earlier  vs  models  of 
Eurasia  (Villasenor  et  a/.,  2001),  but  it  is  not  yet  clear  how  to  assimilate  information  from  Pn 
tomography  into  3-D  models.  This  is  exacerbated  by  the  fact  that  Pn  is  a  very  complicated  phase, 
whose  general  characteristics  exhibit  substantial  regional  variation  (e.g.,  number  and  location  of 
triplications).  It  is,  therefore,  noteworthy  that  the  3-D  vp  models  that  we  construct  from  vs  fit 
regional  P-wave  travel  times  well  and  improve  regional  location  capabilities  relative  to  1-D 
models,  such  as  AK135,  that  have  been  constructed  using  regional  P-wave  data. 

The  quality  of  S- wave  models  of  the  mantle  has  been  improving  rapidly  over  the  past  few  years. 
Continued  improvement,  therefore,  promises  to  advance  location  capabilities  further.  For 
example,  the  model  CUB2.0  was  constructed  by  perturbing  the  Harvard  model  S20A  in  order  to 
fit  a  large  new  data  base  of  regional  group  speed  measurements  in  the  attempt  to  improve  lateral 
and,  more  importantly,  vertical  resolution  in  the  resulting  model.  There  was  no  guarantee  that 
this  procedure  would  produce  a  vp  model  that  would  improve  the  fit  to  regional  P-phases,  but  it 
clearly  does.  Future  advancements  in  vs  models,  both  through  theoretical  developments  as  well  as 
the  assimilation  of  a  variety  of  data  that  have  hitherto  not  been  used  (e.g.,  heat  flow,  Shapiro  & 
Ritzwoller  (2002b);  receiver  functions,  Julia  et  al  (2000),  will  provide  a  natural  path  to  further 
improvements  in  regional  location  capabilities. 

4.3  Importance  and  Promise  of  the  “Theoretical”  vs  to  vp  Conversion 

The  conversion  of  vs  to  vp  has  been  a  major  stumbling  block  to  adopting  global  3-D  seismic 
models  as  part  of  both  teleseismic  and  regional  location  methods.  Recent  advances  in  mineral 
physics  have  greatly  improved  the  theoretical  conversion  of  vs  to  vp.  Continued  improvements  in 
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the  theory  (e.g.,  finite  strain  theory)  as  well  as  developments  in  ancillary  information  (e.g., 
information  about  the  upper  mantle  composition)  can  be  incorporated  naturally  into  this 
conversion  scheme.  In  addition,  the  conversion  can  be  tuned  regionally,  by  varying  mantle 
composition,  and  optimizing  the  anelastic  correction  as  well  as  other  parameters  in  the 
conversion.  To  date,  we  have  applied  the  vs  'to  vp  conversion  directly  from  the  mineral  physics 
literature.  Empirical  phase  path  data  sets,  in  particular,  could  be  used  to  tune  the  conversion  by 
optimizing  certain  parameters  that  are  set  in  the  conversion.  These  factors  argue  that  the 
theoretical  conversion  from  vs  to  vp  provides  a  promising  basis  for  future  improvements  in 
regional  location  capabilities  based  on  3-D  seismic  models. 
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Abbreviations  and  Acronyms 


1-D:  one-dimensional. 

3-D:  three-dimensional. 

AK135:  1-D  model  of  the  Earth  (Kennett  et  al:  1995). 

CTBT:  the  nuclear  Comprehensive  Test  Ban  Treaty. 

CU-Boulder:  University  of  Colorado  at  Boulder. 

CUB2.0:  3-D  vs  model  of  the  Earth  constructed  at  CU-Boulder  (Shapiro  &  Ritzwoller,  2002a). 
CUB2.0JBMP:  3-D  vp  component  of  CUB2.0  model  based  on  empirical  relations  between 
vp  and  vs. 

CUB2.0-TH:  3-D  vp  component  of  CUB2.0  model  based  on  theoretical  relations  between 
vp  and  vs. 

EHB:  Global  earthquake  catalog  (Engdahl,  van  der  Hilst,  and  Buland.  1998). 

DMC:  Data  Management  Center  of  IRIS,  Seattle,  Washington,  U.S. 

GEOFON:  Seismological  Data ‘Center,  the  Federative  Republic  Germany. 

GEOSCOPE:  Seismological  Data  Center,  France. 

GMT:  Generic  Mapping  Tools  (software  package  developed  by  Wessel  &  Smith,  1991,  1995). 

GT:  Ground  Truth  (e.g.  Bondar  et  al ,  2002). 

HDC:  Hypocentroidal  Decomposition  (event  location  method  of  Jordan  Sz  Sverdrup,  1981). 
IASP91:  1-D  Earth  Model  accepted  by  the  International  Association  of  Seismology  and  Physics  of 
the  Earth’s  Interior  (Kennett  &  Engdahl,  1991). 

IDC:  International  Data  Center,  Vienna,  Austria. 

ISC:  International  Seismological  Centre,  United  Kingdom. 

IRIS:  Incorporated  Research  Institutions  for  Seismology. 

NEIS:  U.S.  National  Earthquake  Information  Service,  Golden,  Colorado,  U.S. 

PDE:  Preliminary  Determination  of  Epicenters,  Seismological  Bulletin  issued  by  NEIS. 

PREM:  Preliminary  Reference  Earth  Model  (Dziewonski  Sz  Anderson,  1981). 

S20A:  3-D  vs  mantle  model  constructed  at  Harvard  University  (Ekstrom  &  Dziewonski,  1998). 
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